Dynamical Charge Susceptibility in the Hubbard Model

TL;DR

This paper calculates the dynamical charge susceptibility in the 2D Hubbard model using dynamical cluster approximation, revealing doping dependence and potential for experimental measurement, but finds no clear pseudogap signature.

Contribution

It provides a detailed analysis of charge susceptibility's frequency, doping, and temperature dependence in the Hubbard model, connecting theoretical predictions with experimental observables.

Findings

Susceptibility has a single peak at a characteristic frequency.

Little dependence on momentum or temperature.

Doping influences the static susceptibility significantly.

Abstract

We compute the dynamical charge susceptibility in the two-dimensional Hubbard model within the dynamical cluster approximation. In order to understand the connection between charge susceptibility and pseudogap, we investigate the momentum, doping, and temperature dependence. We find that as a function of frequency, the dynamical charge susceptibility is well represented by a single peak at a characteristic frequency. It shows little momentum or temperature dependence, while the doping dependence is more evident, and no clear signature of the pseudogap is observed. Data for the doping evolution of the static susceptibility and for fluctuation diagnostics are presented. Our susceptibilities should be directly measurable in future Momentum-resolved electron energy-loss spectroscopy experiments.

Equations

$$H=\sum _{k\sigma }({\epsilon }_k-\mu )c_{k\sigma }^{†}{c}_{k\sigma }+U\sum _i{n}_{i\uparrow }{n}_{i\downarrow },$$

$${\chi }_{\sigma {\sigma }^{\prime }}^{\omega {\omega }^{\prime }\Omega }$$

$$\times {\chi }_{\sigma {\sigma }^{\prime }}(k{\tau }_1,(k+q){\tau }_2,(k^{\prime }+q){\tau }_3)d{\tau }_{1}d{\tau }_{2}d{\tau }_3.$$

$${\chi }_d^{\omega {\omega }^{\prime }\Omega }(k,k^{\prime },q)={\chi }_{\uparrow \uparrow }^{\omega {\omega }^{\prime }\Omega }(k,k^{\prime },q)+{\chi }_{\uparrow \downarrow }^{\omega {\omega }^{\prime }\Omega }(k,k^{\prime },q),$$

$${\chi }_d^{\omega {\omega }^{\prime }\Omega }$$

$$\times F_d^{{\omega }_1{\omega }_2\Omega }(k_1,k_2,q){\chi }_0^{{\omega }_2{\omega }^{\prime }\Omega }(k_2,k^{\prime },q),$$

$${\chi }_{ch}(Q,\Omega )=\frac{2}{{\beta }^2{N}^2}\sum _{k,k^{\prime }}\sum _{\omega ,{\omega }^{\prime }}{\chi }_{l,d}^{\omega {\omega }^{\prime }\Omega }(k,k^{\prime },Q).$$

$${\tilde{\Sigma }}_Q(K)$$

$$=$$

Full text

Charge fluctuations in the high-temperature superconductors have generated renewed interest, as new experimental probes, such as Resonant inelastic X-ray scattering RIXS Ishii et al. (2017); Suzuki et al. (2018); da Silva Neto et al. (2018); Hepting et al. (2018); Miao et al. (2019); Ishii et al. (2019) or Momentum-resolved electron energy-loss spectroscopy M-EELS Vig et al. (2017); Husain et al. (2019) may be able to directly measure them as a function of momentum and energy.

Much of the low energy physics of the cuprates is believed to be described by the one-band Hubbard model with an interaction strength close to the bandwidth and a small next-nearest neighbor hopping Anderson (1987); Scalapino (2007). The model has a Mott-insulating state at half filling, as well as pseudogap, superconducting, and Fermi liquid phases in regimes that are remarkably similar to what is seen in the cuprates. Hubbard model simulations of spectral function Huscroft et al. (2001); Civelli et al. (2005); Stanescu and Kotliar (2006); Park et al. (2008); Liebsch and Tong (2009); Ferrero et al. (2009a, b); Sakai et al. (2009, 2010); Lin et al. (2010) show a clear suppression of the density of states at the antinodal (but not at the nodal) points in the pseudogap regime. Results for optical conductivity Ioffe and Millis (2000); Toschi et al. (2005); Millis et al. (2005); Comanac et al. (2008); Ferrero et al. (2010); Lin et al. (2010); Bergeron et al. (2011), Raman spectroscopy Lin et al. (2012); Gull and Millis (2013), neutron spectroscopy LeBlanc et al. (2019) and Nuclear magnetic resonance Macridin et al. (2006); Chen et al. (2017) similarly exhibit the salient features of experiments.

Much less is known about the charge susceptibility of the Hubbard model. Ground state calculations with a variety of methods show various types of charge order at $1/8$ doping White and Scalapino (2003); Hager et al. (2005); Chang and Zhang (2010); Scalapino (2012); Zheng and Chan (2016); Huang et al. (2017); Zheng et al. (2017); Huang et al. (2018); Vanhala and Törmä (2018); Darmawan et al. (2018); Tocchio et al. (2019), competing closely with superconductivity. High-temperature lattice simulations show the doping evolution of the dynamical susceptibility far above pseudogap regime Zimmermann et al. (1997); Sherman and Schreiber (2008); Hochkeppel et al. (2008); Kung et al. (2015). However, no theoretical results are available so far in the pseudogap regime or near the superconducting phase. These results are needed for the attribution of experimental results to pseudogap physics, as experiments often measure a combination of low-energy physics and higher energy contributions from three-dimensional and atomic physics Hepting et al. (2018); Lin et al. (2019).

In this paper, we present a detailed simulation of the momentum, doping, and temperature dependence of the dynamical charge susceptibility $\chi_{ch}(Q,\Omega)$ using eight-site cluster dynamical mean field theory. We also analyze the change of $\chi_{ch}(Q,\Omega)$ in the pseudogap and show the extent to which the establishment of a pseudogap can be understood as a consequence of charge fluctuations.

We solve the two-dimensional Hubbard model with repulsive interaction on a square lattice,

[TABLE]

within the dynamical cluster approximation (DCA) Maier et al. (2005a) on an eight-site cluster. Here $i$ labels the sites in a lattice and $k$ a momentum in the corresponding Brillouin zone. $\varepsilon_{k}=-2t(\cos k_{x}+\cos k_{y})-4t^{\prime}\cos k_{x}\cos k_{y}$ denotes the dispersion, $\mu$ the chemical potential, and $U$ the on-site interaction. We chose $U/t=7$ and $t^{\prime}/t=-0.15$ to differentiate electron from hole doping. Where we present data in K or eV, we use t $\approx$ 0.35 eV corresponding to high-$T_{c}$ superconductors Andersen et al. (1994, 1995); Kim et al. (1998); Hasan et al. (2000). At these parameters, the model is Mott insulating at half filling, and superconducting at low temperature (the maximum $T_{c}$ on the eight-site cluster is near $\mathrm{T}=1/40~{}t$ Gull et al. (2009); Gull and Millis (2012); Gull et al. (2013); Gull and Millis (2014, 2015)), and exhibits a pseudogap region where the single-particle spectral function at the node stays metallic while the antinode is insulating Gull et al. (2009, 2010). Fig. 1 shows an overview of this phase diagram. The DCA can be understood as a momentum-space approximation of the self-energy that coarse-grains the momentum structure, but retains the full frequency dependence Maier et al. (2005a); Fuhrmann et al. (2007). The method is controlled, in the sense that for cluster size $N_{c}\rightarrow\infty$ the exact solution is recovered Maier et al. (2005b); Fuchs et al. (2011); LeBlanc et al. (2015). Our solution is restricted to the paramagnetic phase. Inset of Fig. 1 shows the cluster geometry.

The main result of this paper is the dynamical charge susceptibility $\chi_{ch}(\mathrm{Q},\Omega)$ obtained by solving the DCA equations using a numerically exact continuous-time auxiliary field quantum impurity solver Gull et al. (2008, 2011a, 2011b). The single and two-particle Green’s functions $G_{\sigma}(k_{1}\tau_{1},k_{2}\tau_{2})$ $=$ $\langle T_{\tau}(c_{k_{1}\sigma}^{\dagger}(\tau_{1})c_{k_{2}\sigma}(\tau_{2}))\rangle$ and $G_{2,\sigma_{1}\sigma_{2}\sigma_{3}\sigma_{4}}(k_{1}\tau_{1},...,k_{4}\tau_{4})$ $=$ $\langle T_{\tau}(c_{k_{1}\sigma_{1}}^{\dagger}(\tau_{1})c_{k_{2}\sigma_{2}}(\tau_{2})c_{k_{3}\sigma_{3}}^{\dagger}(\tau_{3})c_{k_{4}\sigma_{4}}(\tau_{4}))\rangle$ define the generalized susceptibility as $\chi_{\sigma\sigma^{\prime}}(k_{1}\tau_{1},k_{2}\tau_{2},k_{3}\tau_{3})$ $=$ $G_{2,\sigma\sigma\sigma^{\prime}\sigma^{\prime}}(k_{1}\tau_{1},k_{2}\tau_{2},k_{3}\tau_{3},0)-G_{\sigma}(k_{1}\tau_{1},k_{2}\tau_{2})G_{\sigma^{\prime}}(k_{3}\tau_{3},0)$, or in Fourier space as Rohringer et al. (2012)

[TABLE]

where $\omega$, $\omega^{\prime}$ are fermionic and $\Omega$ bosonic Matsubara frequencies.

The susceptibility in the density channel is defined as

[TABLE]

and can be decomposed into two parts:

[TABLE]

with $F_{d}$ the full vertex in density channel, $N$ the number of momentum points summed up, and $\chi_{0}$ the bare susceptibility $\chi_{0}^{\omega\omega^{\prime}\Omega}(k,k^{\prime},q)=-\beta NG_{\sigma}(i\omega,k)G_{\sigma}(i\omega+i\Omega,k+q)\delta_{\omega\omega^{\prime}}\delta_{kk^{\prime}}$.

A DCA calculation only yields the “cluster” Green’s functions and susceptibilities $\chi_{c,d}^{\omega\omega^{\prime}\Omega}(K,K^{\prime},Q)$ at the cluster momenta $K,K^{\prime}$ and $Q$. The corresponding lattice susceptibility $\chi_{l,d}(k,k^{\prime},Q)$ is related to this quantity as $\chi_{l,d}^{-1}(k,k^{\prime},Q)-\chi_{l,0}^{-1}(k,k^{\prime},Q)=\chi_{c,d}^{-1}(K,K^{\prime},Q)-\chi_{c,0}^{-1}(K,K^{\prime},Q)$ (for details see Ref. Maier et al., 2005a). The dynamical charge susceptibility is then obtained as

[TABLE]

Eight-site DCA yields results for eight $Q$; four of those are equivalent because of symmetry (see inset of Fig. 1), such that we have four independent momenta: $Q=(0,0)$, $Q=(\pi/2,\pi/2)$, $Q=(\pi,0)$, and $Q=(\pi,\pi)$. We have also performed simulations for 4-site and 16-site clusters at select points and high temperature; these results show qualitatively the same results as 8-site clusters.

We present our simulation data as a function of Matsubara frequency $i\Omega_{n}$. In addition, we also show analytically continued result of $\mathrm{Im}\chi_{ch}(\Omega)$ as a function of real frequency $\Omega$. Analytical continuation Jarrell and Gubernatis (1996); Levy et al. (2017) is an uncontrolled procedure that may exponentially amplify statistical uncertainties, especially in the case of bosonic functions such as the charge susceptibility considered here. Values of $\mathrm{Im}\chi(\Omega)$ at high frequencies are generally less reliable than at low frequency. Our data mostly results in a single large peak at a characteristic frequency, and continuations with different default models do not lead to appreciable differences in this feature. As we will show below, a simple interpretation of this feature in terms of single-particle quantities is not possible due to the importance of vertex functions.

Figure. 2 shows $\chi_{ch}(Q)$ as a function of frequency at T = 0.1t $\sim$ 400 K, at half filling (Figure. 2 a), and in the underdoped (Figure. 2 b), optimally doped (Figure. 2 c), and overdoped (Figure. 2 d) regime. Purple points denote values at $Q=(0,0)$; the value at $i\Omega_{n}=0$ corresponds to the static ($\Omega=0$) uniform ($Q=(0,0)$) charge susceptibility which is small in the insulator and generally rises as doping is increased. A Ward identity Hafermann et al. (2014) requires the frequency dependence to be identically zero in systems that conserve total charge. Unlike many low-order diagrammatic methods Bohm and Pines (1951); Hedin (1965); Bickers (2004); Brener et al. (2008); LeBlanc et al. (2019), DCA satisfies this constraint exactly.

Data for susceptibilities at the three-momentum transfer $(\pi,\pi)$ (blue), $(\pi,0)$ [red, degenerate with $(0,\pi)$], and $(\pi/{2},\pi/{2})$ [orange, degenerate with $(\pm\pi/{2},\pm\pi/{2})$] exhibit a smooth frequency dependence. Remarkably, data at $Q=(\pi/{2},\pi/{2})$ and at $Q=(\pi,0)$ are almost identical, both in the half-filled and in the doped case. The momentum dependence of the charge susceptibility is therefore very different from that of the magnetic susceptibility. For magnetic susceptibility, as found in several approaches LeBlanc et al. (2019); Hochkeppel et al. (2008); Kung et al. (2015), the value at $Q=(\pi,\pi)$ is much larger than any other momentum transfer and rapidly grows as temperature decreases. In contrast, within the momentum resolution achievable within DCA, no dominant contribution to the charge susceptibility is found.

In order to make a connection to experiment, we show analytically continued data corresponding to the Matsubara curves in the lower panel [omitting $Q=(0,0)$]. Within our resolution, our data for $Q=(\pi/{2},\pi/{2})$ and $Q=(\pi,0)$ are described well by a single peak with a maximum near $\omega$ = 7t $\sim$ 2.45 eV (half filling) and 4t $\sim$ 1.40 eV (overdoped). Data at $Q=(\pi,\pi)$ exhibit a peak at a substantially higher frequency (8.5t $\sim$ 2.97 eV for half filling and 5.5t $\sim$ 1.92 eV for overdoped).

Figure. 3 shows the doping dependence of our data at constant temperature T = 0.1t $\sim$ 400 K. The top panels show $Q=(\pi,\pi)$, the middle panels $Q=(0,\pi)$, and the bottom panels $Q=(\pi/2,\pi/2)$. Left panels show Matsubara data, right panels the corresponding analytically continued real frequency data. Four doping points are shown: half filling (purple), underdoped (5% doping, orange), optimally doped (11% doping, red), and overdoped (21% doping, blue).

In Matsubara space, a gradual doping evolution is visible at low frequencies. Zero-frequency values are reduced in comparison to the overdoped values by a factor of about three, while the high-frequency limit remains unchanged. In the real frequency domain, this corresponds to a lowering of the frequency of the charge susceptibility peak, and a general sharpening. As seen previously, no significant momentum dependence is observed, apart from the $\Gamma$ point, which is zero due to charge conservation. This doping evolution is similar to what is found at high temperature Kung et al. (2015).

Figure. 4 shows the temperature dependence of $\chi_{ch}(Q)$ at different doping levels with momentum transfer $Q=(\pm\pi/2,\pm\pi/2)$ at half filling (Figure. 4 a), and in the under-doped (Figure. 4 b), optimally doped (Figure. 4 c), and over-doped (Figure. 4 d) regime. Top panels are the result in Matsubara space, and the bottom panels are the corresponding results in real frequency. Four temperatures are considered here: $\beta$ = 5 (T $\sim$ 800 K, blue), $\beta$ = 7.5 (T $\sim$ 530 K, orange), $\beta$ = 10 (T $\sim$ 400 K, yellow), and $\beta$ = 15 (T $\sim$ 270 K, purple).

As temperature at half filling is decreased, low frequency values are reduced, while the high frequency values remain invariant. In real frequency, the peak of the analytically continued data moves from $\Omega\sim 6t$ to $\Omega\sim 7.5t$. The reduction in low frequency values decreases as we increase doping. All other cases (underdoped, optimally doped, and overdoped) do not show much temperature dependence in both Matsubara and real frequency. The inset of [Figure. 4 (b)] shows the temperature dependence of static, uniform [$i\Omega_{n}=0,Q=(0,0)$] charge susceptibility at doping level $\delta=-0.05$. As shown in Fig. 1, at this doping level, we gradually enter the pseudogap regime as temperature is decreased from $\beta=$ 5 to 15. While Ref. Macridin et al. (2006) shows that $\chi_{m}(Q=(0,0),0)$ is strongly suppressed at these parameters, $\chi_{ch}$ does not show any signature of the pseudogap. The absence of temperature dependence of $\chi_{ch}$ as we enter the pseudogap regime also means that its doping evolution shown in Fig. 3 is not caused by the existence of pseudogap.

It is interesting to compare the values from our calculations to those obtained without vertex corrections, where data are obtained by convolving $\chi_{0}=G*G$. Figure. 5 shows value for the zero-frequency part, as a function of both electron and hole doping. Remarkably, the vertex corrections strongly suppress the overall charge susceptibility and eliminate a large part of the momentum dependence. In particular, the dominant contribution at $(\pi,\pi)$ is reduced to values similar to the other momenta.

An analysis of the frequency dependence (not shown here) shows the biggest discrepancies at the $\Gamma$ point, where the violation of the charge conservation in the absence of vertex corrections leads to a large frequency-dependent contribution for nonzero momenta.

Our data also shows a pronounced dependence of the charge fluctuations on $t^{\prime}$. Whereas the hole-doped side just shows a slowly increasing momentum-independent charge susceptibility and an overall suppression of the $(\pi,\pi)$ susceptibility to the level of the other momenta, the electron-doped side (positive $\delta$, note that approximately $\delta\rightarrow-\delta$ for $t^{\prime}\rightarrow-t^{\prime}$) shows a large enhancement of the (0, 0) susceptibility as compared to other momenta. This region is close to the onset of a first-order coexistence regime Macridin et al. (2006); Gull et al. (2009); Liebsch and Tong (2009) in this model. Notably this is the same parameter region where strong short-range antiferromagnetic fluctuations are present in the magnetic susceptibility.

Two-particle fluctuations such as the charge fluctuations analyzed in this paper are often interpreted as the underlying cause of changes to single-particle observables. From a computational standpoint, both single- and two-particle quantities are computational outcomes of a simulation of Eq. (1). Attributing certain two-particle fluctuations as the “underlying cause” of a change of single-particle features is therefore difficult. However, as defined in Ref. Gunnarsson et al. (2015), it is possible to express the single-particle self-energy, and thereby correlation contributions to the change of the spectral functions, in terms of two-particle quantities, via the (exact) equation of motion. For magnetic and charge fluctuations, these equations are

[TABLE]

with $F_{d}$ defined in Eq. (4), $F_{m}$ the magnetic analog (for detail mathematical form see Ref. [Rohringer et al., 2012]), and $K$, $K^{\prime}$, $Q$ representing pairs of frequency and momentum.

If a single momentum comprises the majority of the self-energy, and contributions from low energies are dominant, then a description in terms of bosonic modes of that type is convenient. This procedure, called “fluctuation diagnostics”, has been successfully applied to non-Fermi liquid Wu et al. (2017); Arzhang et al. (2019) and real-space correlation functions Gunnarsson et al. (2018).

Figure. 6 shows the contributions of fluctuations to the single-particle self-energy at the antinode (left panels) and at the node (right panels) expressed in terms of the charge contributions discussed in this paper (top panels) and in terms of magnetic contributions LeBlanc et al. (2019). As is evident from the lower panels, the pseudogap is well described by short-ranged $Q=(\pi,\pi)$ magnetic fluctuations. A description in terms of charge modes requires similar contributions from all momenta and a much broader frequency range, leading us to conclude that charge fluctuations are not a good way to describe pseudogap physics in the entire parameter range studied here.

In conclusion, we have analyzed the momentum-dependent charge susceptibility in the Hubbard model for a range of dopings and temperatures and for interaction strengths that are thought to be relevant to superconducting cuprates Gull and Millis (2012); Gull et al. (2013). Experimental progress in M-EELS Vig et al. (2017); Husain et al. (2019) and RIXS promises to make this quantity accessible and will provide data that is directly comparable to our results. Our analysis has shown that the dynamical charge susceptibility can be represented by a single peak at a characteristic frequency that exhibits remarkably little momentum dependence, almost no temperature dependence, and a doping dependence that predominantly shifts the peak to lower frequencies as doping is increased. Vertex contributions are essential, as they eliminate the momentum dependence and lower the overall magnitude of the bare susceptibility at $Q=(\pi,\pi)$. As shown by our fluctuation analysis, charge fluctuations are not convenient to understand changes in the single-particle self-energy and spectral function, as terms from all momenta and many frequencies contribute to the self-energy with comparable strength.

While our results are nonperturbative and our solution of the quantum impurity model is numerically exact, the DCA does have important limitations. In particular, due to the limited momentum resolution, DCA is insensitive to stripes with periods larger than our cluster size. Such stripes are seen in experiment Abbamonte et al. (2012); Kohsaka et al. (2007); Fujita et al. (2004); Parker et al. (2010); Lawler et al. (2010); Hashimoto et al. (2010); Daou et al. (2010) and theory White and Scalapino (2003); Hager et al. (2005); Chang and Zhang (2010); Scalapino (2012); Zheng and Chan (2016); Huang et al. (2017); Zheng et al. (2017) in certain parameter regimes, where they compete with the superconducting state, and have been proposed as the cause or a consequence of the pseudogap Lawler et al. (2010); Hashimoto et al. (2010); Daou et al. (2010); Parker et al. (2010).