Doublon-Holon Origin of the Subpeaks at the Hubbard Band Edges

Dynamical mean-field theory (DMFT) studies frequently observe a fine structure in the local spectral function of the SU(2) Fermi-Hubbard model at half filling: In the metallic phase close to the Mott transition, subpeaks emerge at the inner edges of the Hubbard bands. Here we demonstrate that these subpeaks originate from the low-energy effective interaction of doublon-holon pairs, by investigating how the correlation functions of doublon and holon operators contribute to the subpeaks. A mean-field analysis of the low-energy effective Hamiltonian provides results consistent with our DMFT calculation using the numerical renormalization group as an impurity solver. In the SU(3) and SU(4) Hubbard models, the subpeaks become more pronounced due to the increased degeneracy of doublon-holon pair excitations.

Figure 1

Local correlation functions in (a),(b) the metallic and (c) insulating phases of the SU(2) Hubbard model: the local spectral function $A(\omega )$ (blue solid lines), the correlators of doublon $d_{i\nu }$ and holon $h_{i\nu }$ operators [cf. Eq. (1)] (dash-dotted lines), charge susceptibility ${\chi }_c=A_{\delta \widehat{n},\delta \widehat{n}}$ (red dashed lines), and spin (i.e., flavor) susceptibility ${\chi }_s=A_{\overrightarrow{S},\overrightarrow{S}}/3$ (purple dashed lines), with ${\chi }_{c(s)}(\omega )=-{\chi }_{c(s)}(-\omega )$, $A_{d{d}^{†}}(\omega )=A_{h^{†}h}(-\omega )$, and $A_{dh}(\omega )=A_{h^{†}{d}^{†}}(\omega )$. Here $\delta {\widehat{n}}_i\equiv {\widehat{n}}_i-⟨{\widehat{n}}_i⟩$, and ${\overrightarrow{S}}_i$ is the spin operator at site $i$. Each correlator is averaged over different discretizations (see Sec. I B of Ref. [21]), where the corresponding color-matched shaded area provides an estimate for numerical uncertainties, noticeable only in the HBs. Panels (b) and (c) show different solutions for the same value of $U$ in the coexistence regime. In (b), the inset enlarges the region of the QP. We mark the location of spectral features by vertical dotted lines: (a),(b) subpeak position ${\omega }_p$ (defined as the local maximum near the inner HB edge), subpeak width $\delta \omega$ [defined as the minimum positive value satisfying $A({\omega }_p-\delta \omega )=A({\omega }_p)/2$], spin susceptibility peak position ${\omega }_s$, and (c) inner HB edge at $\mathrm{\Delta }/2$, where $\mathrm{\Delta }$ is the Mott gap.

Figure 2

The $U$ dependence of the spectral features: the position ${\omega }_p$ and width $\delta \omega$ of the subpeaks, the peak position ${\omega }_s$ of spin susceptibility ${\chi }_s$, and the Mott gap $\mathrm{\Delta }$ (cf. Fig. 1). Symbols are data points from the $\mathrm{DMFT}+\mathrm{NRG}$ calculations, lines are fits, and shading gives the 95% prediction bounds of fitting. The zeros of the extrapolated fits of $\mathrm{\Delta }$ and ${\omega }_s$ yield estimates for the critical interaction strengths $U_{c1}=2.37(2)$ and $U_{c2}=2.91(1)$, respectively.

Figure 3

(a) Doublon and holon correlators $A_{d{d}^{†}}$ (orange dash-dotted line) and $A_{dh}$ (green dash-dotted line) from our effective theory for the metallic phase. Lower-energy spin dynamics at energies $|\omega |\lesssim {\omega }_s$ and higher-energy scales $|\omega |\gtrsim U/2$ are neglected (as schematically indicated by the gray shading) by employing the generalized SWT together with a mean-field decoupling scheme. $A_{dh}$ is symmetric, while $A_{d{d}^{†}}$ is asymmetric. Both lines have a pair of peaks at $\omega =\pm {\omega }_{dh}$, showing $A_{d{d}^{†}}({\omega }_{dh})>A_{dh}(\pm {\omega }_{dh})>A_{d{d}^{†}}(-{\omega }_{dh})$. This is qualitatively consistent with the $\mathrm{DMFT}+\mathrm{NRG}$ results for $A_{d{d}^{†}}$ and $A_{dh}$ at $\omega =\pm {\omega }_p$ in Fig. 1 using the same color coding. (b) The peak position ${\omega }_{dh}$ from the effective theory decreases linearly with increasing $U$. The narrow shading gives the 95% prediction bounds of a linear fit. ${\omega }_{dh}$ nicely overlaps with ${\omega }_p$ (data taken from Fig. 2) up to an overall scaling factor. We take ${\mathrm{\Delta }}_{dh}=2.91=U_{c2}$ independent of $U$, while the half-filled fraction $⟨P_{i1}⟩$ is $U$ dependent, with the data taken from our $\mathrm{DMFT}+\mathrm{NRG}$ results [21].

Figure 4

Local spectral function $A(\omega )$ for (a) the SU(3) and (b) SU(4) Hubbard models in their metallic phases. Shading again reflects the uncertainties based on discretization averaging (cf. Fig. 1). For $N=3$, the chemical potential $\mu$ was fine-tuned to have the integer filling $⟨{\widehat{n}}_i⟩\simeq 1$ for different $U$, as shown in the legend of (a). For $N=4$, we have $\mu =0$ due to particle-hole symmetry. In all cases, being in the metallic regime, subpeaks emerge at the inner HB edges.