Dynamical effects on superconductivity in BCS-BEC crossover

TL;DR

This study explores how dynamical interactions influence superconductivity during the BCS-BEC crossover, revealing unique features like altered quasiparticle dispersion and phonon softening, with implications for iron-chalcogenide materials.

Contribution

It provides a detailed analysis of dynamical effects on superconductivity in the BCS-BEC crossover using the Holstein model and advanced numerical techniques.

Findings

Absence of dispersion back-bending in Bogoliubov quasiparticles

Non-monotonous dependence of the pairing gap on coupling

Soft phonon spectrum due to Goldstone mode

Abstract

We investigate the dynamical effects of pairing interaction on superconductivity in BCS-BEC crossover by studying the Holstein model at half-filling where the electron-phonon coupling $g$ controls the crossover. The dynamical mean-field theory was employed in combination with the numerical renormalization group technique. The dynamical effects induce distinct features such as absenceof the dispersion back-bending of Bogoliubov quasi-particles, non-monotonous coupling dependence of the pairing gap, and the soft phonon spectrum due to the Goldstone mode of local pair phase fluctuations. Also interesting is the maximum critical temperature being at the normal state phase boundary. Some of these features have intriguing similarities with the recent observations in the FeSe$_{1-x}$S$_x$ and Fe$_{1+y}$Se$_x$Te$_{1-x}$ iron-chalcogenides in the BCS-BEC crossover.

Equations

$$\mathcal{H}=-\frac{t}{\sqrt{q}}\sum _{⟨i,j⟩\sigma }{c}_{i\sigma }^{†}{c}_{j\sigma }+{\omega }_0\sum _i{a}_i^{†}{a}_i+g\sum _i(a_i^{†}+a_i)\left(n_i-1\right),$$

$$U_{eff}(\omega )=\frac{2g^2{\omega }_0}{{\omega }^2-{\omega }_0^2}.$$

$$Re\Delta (\omega ={\Delta }_p)={\Delta }_p,$$

$$T_c={\Lambda }^{-N_{Tc}/2},$$

$$A(ϵ,\omega )=-\frac{1}{\pi }\text{Im}G(ϵ,\omega ),$$

Full text

$Introduction$ – It suits nature of physics to envisage a unified theory to understand the Bardeen-Cooper-Schrieffer (BCS) theory of superconductivity and Bose-Einstein-condensation (BEC) of composite bosons on an equal footing Griffin et al. (1994); Leggett and Zhang (2012); Randeria and Taylor (2014). This BCS-BEC crossover was investigated by Eagles in late 60’s Eagles (1969) and later by Leggett Leggett (1980) and others Nozières and Schmitt-Rink (1985); Randeria et al. (1989); Sá de Melo et al. (1993); Engelbrecht et al. (1997); Chen et al. (2005) by employing the attractive-Hubbard-like models in continuum or discrete lattices. An important conclusion is that the change between the BCS and BEC regimes is gradual so that it is permissible to think of them as two facets of a unified theory. More concretely, the standard picture of the BCS-BEC crossover from the attractive-Hubbard-like models is that the pairing amplitude $\Delta_{p}$ increases monotonically while the superconducting (SC) critical temperature $T_{c}$ increases to have a maximum and decreases as one goes from the BCS to BEC regimes. The change between the regimes can be tuned by varying the ratio $\Delta_{p}/\epsilon_{F}$, where $\epsilon_{F}$ is the Fermi energy, or the strength of the pairing interaction.

Experimental investigations into the crossover have not progressed much because of the difficulties in realizing quantum materials in the BEC regime of $\Delta_{p}/\epsilon_{F}\sim 1$, whereas the ratio is $\sim 10^{-4}$ for conventional superconductors. Then came a burst of research activities because cuprate high-temperature superconductors in the underdoped pseudogap region seem to be in the BEC regime ($\Delta_{p}/\epsilon_{F}\sim 1$) where the carrier density and $\epsilon_{F}$ are strongly suppressed by the Mott physics Lee et al. (2006). However, the cuprate pseudogap seems to exhibit phenomena too much complex to be understood in terms of the BCS-BEC crossover Randeria and Taylor (2014).

This crossover has also been extensively studied in dilute ultracold Fermi alkali gases by tuning the interaction between the Fermi atoms using a Feshbach resonance Gaebler et al. (2010); Regal et al. (2005); Bloch et al. (2008); Chin et al. (2010). This has revealed many insights into the unitary limit which lies in between the BCS and BEC regimes.

Recently, interest in the BCS-BEC crossover in SC materials has been renewed in iron-chalcogenide superconductors because it apparently seems that they can be tuned into both regimes. FeSe doped with isovalent S, FeSe1-xSx, is in the nematic phase for $x\lesssim 0.17$ and in the tetragonal phase for $x\gtrsim 0.17$. The specific heat $C_{p}$ measurements exhibit the classic BCS behavior with a finite jump at $T=T_{c}$ at $x=0$ and 0.1, but, at $x=0.2$, show $C_{p}/T\propto\sqrt{T/T_{c}}$ below $T_{c}$ and the archetypical BEC behavior of the $\lambda$ transition Sato et al. (2018); Mizukami et al. (2018). $T_{c}$ shows a broad maximum and decreases as $x$ changes between 0 and 0.25 Sato et al. (2018); Mizukami et al. (2018). However, the scanning tunneling microscopy measurements show that $\Delta_{p}$ decreases as $x$ is increased into the BEC regime contrary to the standard BCS-BEC crossover picture Hanaguri et al. (2018).

Also interesting is the angle-resolved photoemission spectroscopy study on Fe1+ySexTe1-x by Kanigel group Rinott et al. (2017); Lubashevsky et al. (2012). They altered the Fe content $y$ for a fixed $x=0.4$ to tune the chemical potential, and measured the SC quasi-particle dispersion as the system evolves from BCS to BEC regimes. They found that as $\Delta_{p}/\epsilon_{F}$ increases, the dispersion changes from a characteristic back-bending near $k_{F}$ in BCS regime to a minimum gap at $k=0$ in BEC regime, consistent with the theoretical results of Loh et al Loh et al. (2016). Our study to be presented here indicates this picture needs to be elaborated on.

Now that the crossover is being explored experimentally in the condensed materials and that it is the frequency dependence of the interaction which underlies pairing at least for the conventional superconductors, it seems worthwhile to investigate the frequency dependence of pairing interaction in the BCS-BEC crossover. In this Letter, we build on this idea and show that there are dynamically induced interesting features in the BCS-BEC crossover such as absence of the dispersion back-bending, non-monotonous coupling dependence of the pairing gap, and the maximum $T_{c}$ at the critical coupling concomitantly with the emergence of the Goldstone mode. We will focus on superconductivity in the Holstein model of electron-phonon coupled systems defined in Eq. (1) below to explicitly present these features. They seem to have intriguing similarities with the recent observations in the iron-chalcogenides in the BCS-BEC crossover, yet, our goal is to gain new insights into the BCS-BEC crossover brought in by the dynamical effects rather than to describe specific experiments.

$Model$ – The Holstein model is given by

[TABLE]

where $c_{i\sigma}$ and $a_{i}$ are the fermion of spin $\sigma$ and boson operators at the site $i$, and $\langle i,j\rangle$ implies the nearest neighbors with the coordination number $q$, and $n_{i}=\sum_{\sigma}c_{i\sigma}^{\dagger}c_{i\sigma}$ is the electron density operator at the site $i$. The electrons are coupled with the Einstein phonon of frequency $\omega_{0}$ with the onsite coupling constant $g$. The model has been explored previously via various tools including the dynamical mean-field theory (DMFT) technique Benedetti and Zeyher (1998); Meyer et al. (2002); Capone and Ciuchi (2003); Jeon et al. (2004); Koller et al. (2004).

A useful way to get a grip on this model to integrate the phonons out and write the effective interaction between electrons as

[TABLE]

In the limit $\omega_{0}\rightarrow\infty$, $U_{eff}(\omega)\rightarrow-E_{b}$, where $E_{b}=2g^{2}/\omega_{0}$ is the bipoalron energy. The Holstein model in the regime $\omega_{0}\gg D$ is mapped onto the attractive-Hubbard model with the onsite interaction $|U|=E_{b}$ and the $\omega$ dependence may be neglected. $D=2t$ is the half-bandwidth and is equal to $\epsilon_{F}$ for half-filled cases. Then, the dynamical effects become pronounced for adiabatic cases of $\omega_{0}\ll D$. This regime, however, has been a challenge to treat theoretically because of the different energy scales of electrons and phonons and of emergent superconductivity and soft boson modes Deppeler and Millis (2002). The Holstein model with superconductivity has been calculated previously only for anti-adiabatic regime using DMFT with the continuous-time quantum Monte Carlo technique Murakami et al. (2013, 2014).

Here, we solve the Holstein model at half-filling ($\langle n\rangle=1$) at zero temperature employing DMFT in combination with Wilson’s numerical renormalization group (NRG) technique. The NRG technique can zoom in, via a logarithmic discretization of the conduction band, onto the low energy with an arbitrary precision and is ideal to address the small energy scales in this problem Bauer et al. (2009); Bulla et al. (2008); Jeon et al. (2003); Hewson and Meyer (2002). We will consider $\omega_{0}=0.1$ and $\omega_{0}=2$ as representative cases for adiabatic and anti-adiabatic phonons. We take $D$ as the unit of energy in this paper. Details of the formulation and calculations are given in the Supplementary Material.

*Normal state * – Let us recall what has been established in the normal state Paci et al. (2005); Capone et al. (2006); Yun et al. (2007); Barone et al. (2008). In the absence of long range orders, the model exhibits the 1st order metal-insulator transition as a function of $g$. The insulating phase is a bipolaron insulator (BPI), that is, a local pair state. Like the repulsive or attractive-Hubbard model, there exist two critical values $g_{c1}$ and $g_{c2}$, such that for $g<g_{c1}$ the ground state is a metallic state, for $g>g_{c2}$ a BPI state, and for $g_{c1}\leq g\leq g_{c2}$, both metallic and insulating solutions coexist. We plot in Fig. 1 the $g_{c1}$ and $g_{c2}$ as a function of $\omega_{0}$.

As was reported previously,Yun et al. (2007) $g_{c1}$ and $g_{c2}$ may be estimated by noticing that (a) $U_{c1}/D\approx 2.5$ and $U_{c2}/D\approx 3.0$ for the Hubbard model, and (b) the electron bandwidth is narrowed such that $D\rightarrow De^{-(g/\omega_{0})^{2}}$ because of the polaronic effects Meyer et al. (2002); Benedetti and Zeyher (1998); Capone et al. (2006). We then have $E_{b,c1}/De^{-(g/\omega_{0})^{2}}=2.5$ and $E_{b,c2}/De^{-(g/\omega_{0})^{2}}=3.0$. These estimates for $g_{c1}$ and $g_{c2}$ may be expressed, with the definition $x=2\omega_{0}/U_{c}$ and $y=E_{b,c}/U_{c}$, where the subscript $c$ refers to $c1$ or $c2$, by the single scaling function, $y/\ln(y)=-x$. The numerical results of (a) are plotted on top of this scaling function in the plot (b). Their reasonable overlaps indicate that the polaronic effects play an important role.

*Superconducting $\Delta_{p}$ and $T_{c}$ * – The outcomes of the DMFT calculations are the diagonal self-energy $\Sigma(\omega)$ and off-diagonal self-energy $\phi(\omega)$ in SC state. There is no renormalization of the chemical potential $\mu$ because at the half-filling it is fixed at the center of the band. We will present below the calculations for $T_{c}$, $\Delta_{p}$, density of states, and the spectral function $A(\epsilon,\omega)$ from the obtained self-energies.

The SC gap $\Delta_{p}$ is given by the solution of

[TABLE]

where $\Delta(\omega)=\phi(\omega)/Z(\omega)$ and $Z(\omega)=1-\Sigma(\omega)/\omega$. They are shown in Fig. 2. One can see that for the anti-adiabatic case of $\omega_{0}=2D$ shown in the plot (a), $T_{c}$ increases and have a broad peak and decreases, but $\Delta_{p}$ keeps increasing as $g$ increases like the attractive-Hubbard model. On the other hand, for the adiabatic case of $\omega_{0}=0.1D$ in plot (b), both $T_{c}$ and $\Delta_{p}$ increase and rapidly decrease as $g$ increases. This behavior may be understood by noticing that the bandwidth in the Holstein model is scaled down by the factor of $e^{-(g/\omega_{0})^{2}}$ as discussed above.

$T_{c}$ may be determined by fitting the NRG energy flows ($E_{N}/\Lambda^{-N/2}$) with respect to the NRG site number $N$ as $\Lambda^{(N-N_{Tc})/2}$, and the $T_{c}$ is given by

[TABLE]

where $\Lambda$ is the NRG discretization parameter Hecht et al. (2008). The NRG-DMFT calculations of $T_{c}$ are shown with the black asterisks. See the Supplementary Material for more details. The ratio $2\Delta_{p}/T_{c}$ is also plotted. It shows the BCS value of 3.5 for small $g$ and continues to increase as $g$ increases.

Recall that $T_{c}$ is limited by $\Delta_{p}$ in the BCS and by $D_{S}$ in BEC regimes, where $D_{S}$ is the superfluid stiffness. An estimate of $T_{c}$ taking this observation into account may be done by an interpolation of the gap amplitude and phase stiffness temperature scales of $T_{\Delta}=\Delta_{p}/1.75$ and $T_{D}=D_{S}/2$ such as $T_{c}=(T_{\Delta}^{-\alpha}+T_{D}^{-\alpha})^{-1/\alpha}$. We found that the exponent $\alpha=1.4$ gives the best results from the least square fitting of the above NRG $T_{c}$ calculations. This was represented by the gray dashed curves in Fig. 2 which reasonably reproduce the numerical calculations.

The most interesting for $T_{c}$ is that the maximum occurs at the lower critical value of $g_{c1}$. This holds regardless of the phonon frequency as can be seen from the Fig. 2. It was already noticed in the attractive Hubbard model by Toschi et al. that the maximum $T_{c}$ occurs at $U_{c1}$ Toschi et al. (2005); Kuleeva et al. (2014). One may expect that $T_{c}$ becomes maximum around the crossover regime where the $T_{c}$ limiting temperature scales $T_{\Delta}$ and $T_{D}$ become comparable. The present calculations show that the maximum $T_{c}$ occurs precisely at $g_{c1}$ where appears the soft phonon mode associated with the bipolaron instability.

As the electron-phonon coupling $g$ increases the phonon frequency is renormalized and a zero frequency mode appears at the critical coupling. As $g$ is further increased it hardens back in normal state Jeon et al. (2004). In SC state, however, for $g$ larger than the critical coupling (BEC regime) the zero frequency weight remains non-vanishing. This is due to the Goldstone mode of the phase fluctuations of local pairs Micnas et al. (1990). The last term in the Holstein model of Eq. (1) describes the coupling between the phonons and the local electron density. The phonon spectrum therefore reflects the local pairs in the BEC state. Fuller discussions on the Goldstone mode, phonon spectra, and their connection with the maximum $T_{c}$ will be reported separately.

$Dynamical~{}properties$ – The electronic density of states in SC states are presented in Fig. 3. Plot (a) is for the anti-adiabatic and (b) is for the adiabatic cases. If the system is metallic when SC is turned off (BCS regime), then DOS exhibits only one peak (for $\omega>0$) at $\omega=\Delta_{p}$ in SC state. On the other hand, if insulating when SC is turned off (BEC regime), DOS may exhibit two peaks in SC state: one at $\omega=\Delta_{p}$ from SC and the other around $\omega\approx E_{b}/2$ from the insulating gap. These two peaks clearly show up only around crossover regime for the adiabatic case as can be seen from the plot (b) for $g/D=0.2$. For other regimes only single peak shows up in DOS either because the two energy scales merge into one (anti-adiabatic case) or the pairing coherence peak is suppressed due to the small superfluid stiffness (adiabatic case).

The phase coherence is due to the superfluid stiffness, and we already saw that $D_{S}$ is significantly suppressed in the BEC regime of the adiabatic case from Fig. 2(b). However, $D_{S}$ does not decrease rapidly in the anti-adiabatic case as shown in Fig. 2(a). Therefore, as shown in Fig. 3(a), the coherence peaks remain unsuppressed together with the local pair states at $\omega\approx\Delta_{p}\approx E_{b}/2$ in BEC regime. The higher order peaks also appear at the addition of phonon energy $\omega_{0}$.

$Spectral~{}intensity~{}A(\epsilon,\omega)$ – We now show the spectral function in Fig. 4.

[TABLE]

where, $\epsilon$ stands for $\epsilon_{k}$. The spectral intensity plot in the plane of of $\epsilon$ and $\omega$ represents the quasiparticle dispersion relation, which is probed by the angle-resolved photoemision spectroscopy (ARPES) experiments. We present the representative spectral intensity plots for BCS, crossover, and BEC regimes, respectively, in Fig. 4(a), (b), and (c) for the adiabatic case, and in (d), (e), and (f) for the anti-adiabatic case. The bare dispersion is along the line of $\omega=\epsilon$. In the BCS regime, the back-bending behaviors around the Fermi surface ($\epsilon=0$) can clearly be seen for both adiabatic (plot (a)) and anti-adiabatic (plot (d)) cases. The kink around the bare phonon frequency ($\omega\sim\omega_{0}$) is also observable for the adiabatic phonon in accord with the mean-field theory calculation of Sandvik et al Sandvik et al. (2004).

In the BEC regime, the quasi-particle dispersion for adiabatic and anti-adiabatic cases behave differently as shown in Fig. 4(c) and (f). For the adiabatic case, the dispersion curve, as $\epsilon$ increases, passes through the Fermi surface without the back-bending and resumes the normal state dispersion ($\omega=\epsilon$). The presumed back-bending feature (in, for example, $\epsilon>0$ and $\omega<0$ quadrant) is absent which appears in all regimes in anti-adiabatic case. Instead, the dispersion gap keeps decreasing as one moves away from the bare Fermi surface. The absence of the back-bending was reported by Rinott et al by the ARPES experiments on $\text{Fe}_{1+y}\text{Se}_{x}\text{Te}_{1-x}$ as mentioned in introduction Rinott et al. (2017).

$Summary~{}and~{}discussion$ – Our main results for dynamical effects on the crossover are: (a) The pairing gap $\Delta_{p}$ does not keep increasing; $\Delta_{p}$ increases and decreases as the coupling $g$ increases for adiabatic phonon cases. (b) The maximum $T_{c}$ occurs at the lower critical coupling $g_{c1}$ concomitantly with the emergence of the Goldstone mode due to the local pair phase fluctuations. (c) The SC quasi-particle dispersion in the crossover or BEC regime does not exhibit the back-bending in the adiabatic case. The (b) holds regardless of the phonon frequency, and (a) and (c) are pronounced in the adiabatic cases.

It is important to note that the quasi-particle dispersion may or may not exhibit the back-bending in the BEC regime of superconductivity from dynamically induced insulators. In contrast, the back-bending always shows up in superconductivity emerging from insulators induced by an instantaneous interaction like the attractive-Hubbard model, while it is absent in superconductivity from band insulators Loh et al. (2016). The absence of SC back-bending in latter is understood in the simple mean-field theory picture as follows: Suppose a band insulator where the chemical potential is below the bottom of a conduction band and the minimum gap position is at $k=0$. In SC state, this is in BEC regime and the dispersion gap minimum also occurs at $k=0$. The topology of the minimum gap loci is a finite volume contour for BCS but a point in ${\bf k}$-space in BEC regimes. However, in the present problem the chemical potential is always at the center of a particle-hole symmetric band to ensure the half-filling, and above picture of the chemical potential out of a band can not be applied.

Superconductivity emerging from a dynamically induced insulator we consider here, the absence of back-bending can be traced back to the insulating behavior of the diagonal self-energy $\Sigma(\omega)$ in SC state, that is, the slope of $Re\Sigma(\omega)\rightarrow\infty$ as $\omega\rightarrow 0$. In both adiabatic and anti-adiabatic cases the self-energy in the local pair insulating regime shows this behavior in the normal state like in the Mott insulating state Bulla et al. (2001). In SC state, however, the low energy fluctuations leading to insulating behavior are blocked out in anti-adiabatic cases by the large energy gap of $\Delta_{p}$. But this blocking is not complete in adiabatic cases and the insulating behavior survives. Like the band insulator cases, this insulating behavior underlies the absence of the back-bending. See the Supplementary Material for more details. Application of this idea to Fe1+ySexTe1-x needs more considerations. In particular it seems that the multigap nature of the material ought to be considered in the context of dynamical effects on superconductivity. It will be interesting to observe the BCS-BEC crossover out of interaction induced insulators.

$Acknowledgment$ – We acknowledge the supports from Samsung Science & Technology Foundation through SSTF-BA1502-06, and from National Research Foundation of Korea through NRF-2016R1A6A3A11930653 (THP) and NRF-2018R1D1A1B07043997 (HYC).