Two-particle correlation functions in cluster perturbation theory: Hubbard spin susceptibilities

Cluster perturbation theory (CPT) is a computationally economic method commonly used to estimate the momentum- and energy-resolved single-particle Green's function. It has been used extensively in direct comparisons with experiments that effectively measure the single-particle Green's function, e.g., angle-resolved photoemission spectroscopy. However, many experimental observables are given by two-particle correlation functions. CPT can be extended to compute two-particle correlation functions by approximately solving the Bethe-Salpeter equation. We implement this method and focus on the transverse spin susceptibility, measurable via inelastic neutron scattering or with optical probes of atomic gases in optical lattices. We benchmark the method with the one-dimensional Fermi-Hubbard model by comparing with known results.

Figure 1

Diagrams used in cluster perturbation theory: (a) the Dyson equation, (b) the Bethe-Salpeter equation for the generalized four-point susceptibility, (c) the Bethe-Salpeter equation for the two-point susceptibility, and (d) the analog of the Dyson equation for the two-point susceptibility. The four-point vertex, $\mathrm{\Gamma }(i\nu ,i{\nu }^{\prime },i\omega )$, is approximated by the two-point vertex, $\mathrm{\Gamma }\left(i\omega \right)$. Dashed lines represent the noninteracting Green's function $G_0$, while double lines represent the interacting Green's function $G$.

Figure 2

CPT approximation [Eq. (3) to the Green's function, $-ImG(k,\nu +i\eta )$, for the one-dimensional Hubbard model at half filling plotted as a function of wave vector $k$. We choose a cluster size $L=16$, small broadening parameter $\eta =0.2$, and plot $N_p=144$ points for $k$ and $\nu$. From top to bottom, we have Hubbard interaction strengths: $U=0,2,4,$ and 8.

Figure 3

Left: CPT approximation [Eq. (10) to the spin susceptibility, $Im\chi (q,\omega +i\eta )$ for the one-dimensional Hubbard model at half filling plotted as a function of wave vector. We choose $L=16,\eta =0.5$, and $N_p=96$. Right: The same but for the RPA-CPT approximation, Eq. (11). The interaction is chosen to be weak, $U=0.1$ and $U=0.5$, (top and bottom, respectively) so the RPA-CPT is accurate. The oscillating lines are a numerical artifact discussed in Appendix pp5.

Figure 4

Same as Fig. 3 but for an intermediate interaction strength, $U=1$ (top) and $U=2$ (bottom). The CPT (left) and RPA-CPT (right) agree well up to $U=2$ where deviations between the methods start to appear.

Figure 5

Left: The same as Fig. 3 but for larger interaction strengths, $U=4$ and $U=8$, and with $\eta =0.2$. Right (top): DMRG result for a 64 site chain and $U=4$. Here we see that the CPT method broadens the peak at $q=\pi$ in comparison to DMRG. Right (bottom): Müller estimate, Eq. (12), for $U=8$ with a cutoff intensity such that $\text{max}\left[{\chi }_{\text{Müller}}(q,\omega )\right]=\text{max}[Im{\chi }_{\text{CPT}}(q,\omega +i\eta )]$.

Figure 6

Left: CPT approximation [Eq. (10) to the spin susceptibility, $Im\chi (q,\omega +i\eta )$ for the one-dimensional Hubbard model at one-quarter filling (density $n=1/2)$ plotted as a function of wavevector. We choose $L=16,\eta =0.2$, and $N_p=96$. The chemical potentials used are listed in Appendix pp8. Numerical artifacts in CPT are due to pole mismatches discussed in Appendices pp4 and pp5. Right: DMRG result for a 64-site chain at the same filling with $\eta =.05$.

Figure 7

Same as Fig. 6 except for one-eighth filling (density $n=1/4)$.

Figure 8

Top: CPT susceptibility for $L=16,U=0.5$ with $\eta =0.2$ to highlight the numerical instabilities arising from a low broadening parameter. The white dots appear where $\omega <0$. The badly behaved points follow along an oscillating structure with peaks proportional to cluster size $L$. Bottom: The same but for $U=1.0$. The oscillating structure shrinks as the misbehaving poles lose weight at lower energy scales.

Figure 9

Frequency dependence of the CPT two-point vertex function $\mathrm{\Gamma }(q,\omega )$ used in the Bethe-Salpeter equations [see Fig. 1 for $H_{\text{1D}}$. The horizontal dashed lines show the RPA approximation $\mathrm{\Gamma }(q,\omega )=U$. The solid lines show the CPT approximation obtained from inverting Eqs. (7a) with $\eta =0.5$ and $N_p=48$. Here we see that at low $U$, the real part of the CPT approximation at low frequencies is very close to the RPA. The imaginary part of the CPT vertex function is also shown for comparison.

Figure 10

Finite-size trend of the gap $\mathrm{\Delta }$. $\mathrm{\Delta }$ is obtained as the frequency for which $\chi (\pi ,\omega )$ is first nonzero, using CPT for $U=4$ (top) and $U=8$ (bottom) in $H_{\text{1D}}$. Here we see a trend toward zero. Linear fits of the last four data points show $y$ intercepts of 0.0255(8) and 0.020(1) for $U=4$ and $U=8$, respectively.