Quasiparticle Interference and Symmetry of Superconducting Order Parameter in Strongly Electron-Doped Iron-based Superconductors

TL;DR

This paper investigates the pairing symmetries in electron-doped iron-based superconductors with strong spin-orbit coupling, identifying specific states consistent with experimental quasiparticle interference data and proposing their microscopic realization.

Contribution

It analyzes quasiparticle interference signatures to distinguish pairing symmetries, highlighting the relevance of spin-triplet dominated s-wave and d-wave states in these materials.

Findings

Only certain pairing states match experimental QPI data.

Spin-orbit coupling leads to mixed singlet-triplet pairing components.

Proposed microscopic models support the realization of these pairing states.

Abstract

Motivated by recent experimental reports of significant spin-orbit coupling (SOC) and a sign-changing order-parameter in the Li$_{1-x}$Fe$_x$(OHFe)$_{1-y}$Zn$_y$Se superconductor with only electron pockets present, we study the possible Cooper-pairing symmetries and their quasiparticle interference (QPI) signatures. We find that each of the resulting states - $s$-wave, $d$-wave and helical $p$-wave - can have a fully gapped density of states (DOS) consistent with angle-resolved photoemission spectroscopy (ARPES) experiments and, due to spin-orbit coupling, are a mixture of spin singlet and triplet components leading to intra- and inter-band features in the QPI signal. Analyzing predicted QPI patterns we find that only the spin-triplet dominated even parity $A_{1g}$ (s-wave) and $B_{2g}$ (d-wave) pairing states are consistent with the experimental data. Additionally, we show that these states can indeed be realized in a microscopic model with atomic-like interactions and study their possible signatures in spin-resolved STM experiments.

Tables

Table 1. Table 1: Fitting parameters relevant for intercalated FeSe taken from Ref. Eugenio and Vafek, 2018 .

${\alpha }_1$			$782.512\text{meV}{\mathring{A}}^2$
${\alpha }_3$			$-1400\text{meV}{\mathring{A}}^2$
$\frac{1}{2m_1}$			$-492.01\text{meV}{\mathring{A}}^2$
$\frac{1}{2m_3}$			$1494.14\text{meV}{\mathring{A}}^2$
$v$			$224.406\text{meV}\mathring{A}$
${\lambda }_z$			$26\text{meV}$
$p_{z_1}=p_{z_3}$			$0$

Equations

$$H=\sum _{\sigma ,{\sigma }^{\prime }=\uparrow ,\downarrow }\sum _{\text{k}}$$

$$\displaystyle\Psi_{X,\sigma}(\textbf{k})=\left(\begin{array}[]{c}c_{yz,\sigma}(\textbf{k})\\ c_{xy^{X},\sigma}(\textbf{k})\end{array}\right)\quad,\quad\Psi_{Y,\sigma}(\textbf{k})=\left(\begin{array}[]{c}c_{xz,\sigma}(\textbf{k})\\ c_{xy^{Y},\sigma}(\textbf{k})\end{array}\right)$$

$$h_X^{\prime \sigma {\sigma }^{\prime }}(\mathbf{k})$$

$$\displaystyle+p_{z_{1}}(k_{x}^{3}-k_{y}^{3})+p_{z_{2}}k_{x}k_{y}(-k_{x}+k_{y})\Big{)}\sigma^{z}_{\sigma,\sigma^{\prime}}\tau_{1}$$

$$h_Y^{\prime \sigma {\sigma }^{\prime }}(\mathbf{k})$$

$$\displaystyle+p_{z_{1}}(k_{x}^{3}+k_{y}^{3})+p_{z_{2}}k_{x}k_{y}(k_{x}+k_{y})\Big{)}\sigma^{z}_{\sigma,\sigma^{\prime}}\tau_{1}$$

$$h_X(\mathbf{k})$$

$$h_Y(\mathbf{k})$$

$${\Lambda }_{\text{SOC}}^{\sigma ,{\sigma }^{\prime }}=i{\lambda }_{\text{SOC}}\left(\frac{{\tau }_1+i{\tau }_2}{2}\otimes {\sigma }_{\sigma ,{\sigma }^{\prime }}^x+\frac{{\tau }_1-i{\tau }_2}{2}\otimes {\sigma }_{\sigma ,{\sigma }^{\prime }}^y\right).$$

$$H_{\text{int}}(\text{r})=\sum _{\alpha }$$

$$+$$

$$d_{xz,\sigma }^{A(B)}(\text{k})$$

$$d_{yz,\sigma }^{A(B)}(\text{k})$$

$$d_{xy,\sigma }^{A(B)}(\text{k})$$

$$H_{\text{int}}^{A_{1g}}+H_{\text{int}}^{B_{2g}}+H_{\text{int}}^{E_g}$$

$$\displaystyle=\sum_{\textbf{k},\textbf{q}}\Big{(}\Psi^{\dagger}_{1\uparrow}(\textbf{k})\tau_{0}\Psi^{*}_{1\downarrow}(-\textbf{k}),\Psi^{\dagger}_{3\uparrow}(\textbf{k})\tau_{0}\Psi^{*}_{3\downarrow}(-\textbf{k})\Big{)}\times$$

$$\displaystyle\times\left(\begin{array}[]{cc}\frac{1}{2}(U+J_{11}^{\prime})&J_{13}^{\prime}\\ J_{13}^{\prime}&U\\ \end{array}\right)\left(\begin{array}[]{c}\Psi^{T}_{1\downarrow}(-\textbf{q})\tau_{0}\Psi_{1\uparrow}(\textbf{q})\\ \Psi^{T}_{3\downarrow}(-\textbf{q})\tau_{0}\Psi_{3\uparrow}(\textbf{q})\end{array}\right)$$

$$+$$

$$A_{1g}^s:$$

$${\Delta }_3^A{\Psi }_{3,\sigma }^T(-\text{k}){\tau }_{0}i{\sigma }^y{\Psi }_{3,{\sigma }^{\prime }}(\text{k})$$

$$B_{2g}^s:$$

$${\Delta }_3^B{\Psi }_{3,\sigma }^T(-\text{k}){\tau }_{3}i{\sigma }^y{\Psi }_{3,{\sigma }^{\prime }}(\text{k}).$$

$$\displaystyle E_{g1}^{-}=i\left(\begin{array}[]{cc}&\tau_{-}\\ -\tau_{+}&\end{array}\right)\quad E_{g2}^{-}=i\left(\begin{array}[]{cc}&\tau_{+}\\ -\tau_{-}&\end{array}\right)$$

$$A_{1g}^t:{\Delta }_t^A{\Psi }_{M,\sigma }^T(-\text{k})(-E_{g1}^{-},i{E}_{g2}^{-})({\sigma }^z,{\sigma }_o){\Psi }_{M,\sigma }(\text{k})$$

$$B_{2g}^t:{\Delta }_t^B{\Psi }_{M,\sigma }^T(-\text{k})(-E_{g1}^{-},-i{E}_{g2}^{-})({\sigma }^z,{\sigma }_o){\Psi }_{M,\sigma }(\text{k})$$

$$\displaystyle H_{BdG}^{A_{1g}(B_{2g})}=\sum_{\textbf{k}}\Psi^{\dagger}_{g}(\textbf{k})\left(\begin{array}[]{cc}\mathcal{H}_{0}(\textbf{k})&\hat{\Delta}_{A_{1g}(B_{2g})}\\ \hat{\Delta}^{\dagger}_{A_{1g}(B_{2g})}&-\mathcal{H}_{0}(\textbf{k})\end{array}\right)\Psi_{g}(\textbf{k})$$

$$\displaystyle\mathcal{H}_{0}(\textbf{k})=\left(\begin{array}[]{cc}h^{\prime\uparrow\uparrow}_{X}(\textbf{k})&\lambda_{\text{SOC}}\Lambda\\ \lambda_{\text{SOC}}\Lambda^{\dagger}&h^{\prime\downarrow\downarrow}_{Y}(\textbf{k})\end{array}\right)\quad,\quad\Lambda=\left(\begin{array}[]{cc}0&i\\ 1&0\end{array}\right)$$

$${\hat{\Delta }}_{A_{1g}}$$

$${\hat{\Delta }}_{B_{2g}}$$

$${\Psi }_{\text{k}}^{†}=\left({\Psi }_{M,\uparrow }^{†}(\text{k}),{\Psi }_{M,\downarrow }^{†}(\text{k}),{\Psi }_{M,\downarrow }^T(-\text{k}),-{\Psi }_{M,\uparrow }^T(-\text{k})\right).$$

$$\displaystyle H=\frac{1}{2}\sum_{\mathbf{k}}\Psi^{\dagger}(\mathbf{k})\left(\begin{array}[]{cc}\mathcal{H(\mathbf{k})}&\hat{\Delta}\\ \hat{\Delta}^{\dagger}&-\mathcal{H(\mathbf{k})}\end{array}\right)\Psi(\mathbf{k}),$$

$$\displaystyle\mathcal{H(\mathbf{k})}=\left(\begin{array}[]{cccc}h^{\prime\uparrow\uparrow}_{X}(\textbf{k})&&&\Lambda\\ &h^{\prime\uparrow\uparrow}_{Y}(\textbf{k})&-\Lambda^{T}&\\ &-\Lambda*&h^{\prime\downarrow\downarrow}_{X}(\textbf{k})&\\ \Lambda^{\dagger}&&&h^{\prime\downarrow\downarrow}_{Y}(\textbf{k})\end{array}\right).$$

$$\rho (\omega )$$

$$H_{\text{imp}}^{1Fe}(\text{k},{\text{k}}^{\prime })=V_{\alpha ,\beta }{d}_{\alpha ,\sigma }^{A†}(\text{k})d_{\beta ,\sigma }^A({\text{k}}^{\prime }).$$

$$H_{\text{imp}}^{2Fe}(\text{k},\text{k}\prime )={\Psi }_{M,\sigma }^{†}\tilde{V}{\Psi }_{M,\sigma }$$

$$\displaystyle\tilde{V}=\left(\begin{array}[]{cccc}V_{yz,yz}&&&\\ &V_{xy^{X},xy^{X}}&&V_{xy^{X},xy^{Y}}\\ &&V_{xz,xz}&\\ &V_{xy^{X},xy^{Y}}&&V_{xy^{Y},xy^{Y}}\end{array}\right).$$

Full text

I Introduction

In iron-based superconductors, it has been widely believed that superconductivity is driven by repulsive interactions, enhanced by the presence of the spin fluctuations associated with the parent antiferromagnetic state. In this scenario, these fluctuations drive a sign reversal ($s_{\pm}$ state) Kuroki et al. (2008); Mazin et al. (2008); Chubukov et al. (2008); Chubukov (2012); Hirschfeld (2016) between order parameters (OP) on the electron and hole Fermi surface pockets at the $M$- and $\Gamma$-point, respectively. The discovery of superconductivity in intercalated or monolayer FeSe at a critical temperature of the order above 40K revived interest in Fe-based superconductivity, but raised further questions on the origin of superconductivity in these compoundsGuo et al. (2010); Lu et al. (2015); Wang et al. (2012); Miyata et al. (2015); Zhao et al. (2016), because, unlike bulk FeSe, ARPES experiments show that many of these FeSe-derived systems appear to be missing the hole pockets at the $\Gamma$-point required in the the conventional scenario.

Initial model calculations based on the multiorbital spin-fluctuation framework for systems manifesting only electron pockets at the $M$-point predicted $d$-wave symmetry state in this caseMaier et al. (2011); Wang et al. (2011), driven by the spin fluctuations connecting the electron pockets that remain when the hole pockets are removed. In the proper 2-Fe unit cell, such a state must have gap nodes on the Fermi surfaceMazin (2011). This is because the electron pockets located near $(\pi,0)$ or $(0,\pi)$ points of the Brillouin Zone (BZ) in the 1-Fe unit cell fold onto $(\pi,\pi)$ point of the folded BZ as the crystallographic symmetry lowers due to the Se positions. This may lead to hybridization between the electron pocketsKhodas and Chubukov (2012), which then forces the $d_{x^{2}-y^{2}}$-state to acquire gap nodes, although in principle the nodal area may be very small, proportional to the hybridization (“quasinodes”). On the other hand, ARPES experiments in most of the electron-intercalated materials indicated a nodeless superconducting (sc) stateYan et al. (2016); Zhang et al. (2016). Several proposals for the gap structure have been put forward, including a conventional $s^{++}$-wave scenario based on the electron-phonon interaction and orbital fluctuationsOnari and Kontani (2012), as well as the “bonding-antibonding” scenarioMazin (2011); Hirschfeld et al. (2011) in which the order parameter on the inner electron pocket (mostly $d_{xz}/d_{yz}$ character) has one sign, and on the outer electron pocket (mostly $d_{xy}$ character) the otherKhodas and Chubukov (2012). Furthermore, it has been argued that the hybridization of the electron pockets is mainly due to SOCKreisel et al. (2013); Cvetkovic and Vafek (2013); Eugenio and Vafek (2018), which within a 3D spin fluctuation framework may stabilize the bonding-antibonding $s_{\pm}$ state against d-waveKreisel et al. (2013) and mixes a spin-triplet component into the even parity $s^{+-}$-wave stateCvetkovic and Vafek (2013); Eugenio and Vafek (2018). Overall, one can see that the sign structure of the superconducting order parameter is intimately related to the pairing mechanism. Therefore, experiments allowing to determine it could be of great potential importance.

One rapidly developing technique to determine the phase structure of the order parameter makes use of QPI as measured by Fourier transform scanning tunneling microscopy (FT-STM). This probe measures the wavelengths of Friedel oscillations caused by impurities present in a metallic or superconducting system, which in turn contains information on the electronic structure of the pure system. A subset of scattering wave vectors q corresponding to peaks in the FT-STM can be enhanced or not according to the type of disorder and the phase structure of the superconducting gapPereg-Barnea and Franz (2003); Nunner et al. (2006). Recently it was proposed by Hirschfeld, Altenfeld, Eremin and Mazin (HAEM)Hirschfeld et al. (2015) that the sign structure of the order parameter in a multiband system can be extracted from the Fourier transform of the local density of states QPI pattern near an impurity in the superconducting state. The antisymmetrized QPI intensity integrated over the wavevectors corresponding to scattering between two bands was shown to have a dependence on frequency very different for sign-changing and sign preserving scenarios leading to a strong, single-sign enhancement of the integrated response in the former case. This qualitative result was also confirmed by extensive numerical simulations with finite disorderMartiny et al. (2017). Recently, a complementary phase sensitive technique to detect sign-changing gaps in the presence of strong impurity bound states was proposedChi et al. (2017).

Using QPI analysis, the authors of Ref. Sprau et al., 2017 were able to identify a sign changing order parameter in FeSe. Most importantly for our purposes, similar conclusions were recently drawn for the strongly electron doped iron-based superconductor lithium hydroxide intercalated FeSeDu et al. (2017). In other words, the order parameter in Li1-xFex(OHFe)1-yZnySe, alternates sign, either between the Fermi surface sheets, or within individual sheets. However, distinguishing between these alternatives was beyond the resolution of the experiment. In any case, the situation is somewhat more complicated than anticipated in Ref. Du et al., 2017, since the effect of spin-orbit interaction on pairing needs to be taken into account as well. Moreover, recent observation of Majorana zero modes in (Li0.84Fe0.16)OHFeSeLiu et al. (2018) suggests possible broader implications of the spin-orbit coupling for the Cooper-pairing in electron doped intercalated iron-based superconductors. Note that in contrast to Ref. Du et al., 2017, no Zn substitution was used in Ref. Liu et al., 2018. The amount of Zn, however, is relatively small (less than 2 percent). This amount does not affect the superconducting transition temperature or the electronic structure in a significant way and is done only for the purpose of enhancing the QPI signal in the scanning tunneling microscopy.

In this manuscript we study the possible Cooper-pairing symmetries and their QPI signatures for strongly electron-doped Fe-based superconductors using the effective three-orbital model of Refs. Cvetkovic and Vafek, 2013; Eugenio and Vafek, 2018 with spin-orbit coupling and proper consideration of all lattice symmetries of the FeSe space group. We find that each of the resulting states - $A_{1g}$-wave, $B_{2g}$-wave and helical $E_{u}$-wave - can have a fully gapped DOS consistent with ARPES experiments and, due to spin-orbit coupling, are a mixture of spin singlet and triplet components leading to intra- and inter-band features in the QPI signal. Analyzing predicted QPI patterns we find that $A_{1g}$-wave pairing state, with the two dominant peaks in the DOS roughly corresponding to the gap energies on each pocket, and $B_{2g}$-wave pairing state both with a significant even parity spin triplet component are consistent with the experimental data. Moreover, we show that pairing states with dominant spin triplet component can be identified using spin-resolved STM.

II Model

We wish to describe the low energy states near the M,-points of the Brillouin zone using the orbitally projected band model of Ref. Cvetkovic and Vafek, 2013; Eugenio and Vafek, 2018 for the two-iron unit cell. Near the Fermi level, only the $xz$, $yz$ and $xy$ orbitals contribute significantly; hence, the full 10-orbital tight binding model is projected onto the subspace of these three orbitals. The effective low energy Hamiltonian near the M-point that takes into account all the lattice symmetries of the FeSe space group as well as time reversal symmetry is defined as

[TABLE]

where the four component spinor $\Psi^{\dagger}_{M,\sigma}(\textbf{k})=\left(\Psi^{\dagger}_{X,\sigma}(\textbf{k}),\Psi^{\dagger}_{Y,\sigma}(\textbf{k})\right)$ describes the states at the M-point for each spin projection $\sigma$. The doublets $\Psi_{X,\sigma}(\textbf{k})$ and $\Psi_{Y,\sigma}(\textbf{k})$ are defined as

[TABLE]

Moreover, we have

[TABLE]

and

[TABLE]

where the Pauli matrices $\{\sigma^{x},\sigma^{y},\sigma^{z}\}$ and $\{\tau_{1},\tau_{2},\tau_{3}\}$ act on spin and orbital space, respectively. The $\lambda_{z}$, $p_{z_{1}}$ and $p_{z_{2}}$ terms in eq.(LABEL:InterBandSoC) describe the k-dependent intra-band SOC which does not couple the two Fermi pockets but lifts the out-of plane spin degeneracy. The inter-band SOC term which hybridizes X- and Y-pocket is given by

[TABLE]

In order to describe intercalated FeSe we use the Luttinger invariants, Tab.(1), which were evaluated in Ref.Eugenio and Vafek, 2018 based on the available ARPES data.

The value of the intraband spin-orbit coupling, $\lambda_{z}$ is in agreement with those found in ab-initio calculations and ARPES experimentsBorisenko et al. (2016); Day et al. (2018). Furthermore, the value of the interband spin-orbit coupling, $\lambda_{\text{SOC}}$ between the electron pockets separated by the large momentum yielding their hybridization and splitting on the Fermi surface is found to be smaller and is taken to be $\sim 5$ meVDay et al. (2018). Diagonalizing Eq.(3) yields four bands: two regular ones that form the inner and outer electron pocket, see Fig.1(b), and two incipient bands that do not cross the Fermi level. The effect of inter-band SOC on the band structure is visualized in Fig.1(c) and Fig.1(d).

III Mean field phase diagram

Although phenomenologically the classification of superconducting orders for two electron pockets was considered previouslyCvetkovic and Vafek (2013); Eugenio and Vafek (2018) we analyze here its microscopic formulation via mean-field treatment of the atomic on-site interactions given by the Hubbard and Hund’s couplings $U$, $U^{\prime}$, $J$ and $J^{\prime}$ which enter the Hubbard-Kanamori Hamiltonian as

[TABLE]

Here $\alpha\in\{A,B\}$, $\{\sigma,\sigma^{\prime}\}\in\{\uparrow,\downarrow\}$ and $\{\mu,\nu\}\in\{yz,xz,xy\}$ label lattice sites, spins and orbitals, respectively. $d^{\alpha\dagger}_{\mu,\sigma}(\textbf{r})$ and $d^{\alpha}_{\mu,\sigma}(\textbf{r})$ are the second quantized operators creating and annihilating particles on sub-lattice A and B, see Fig.1(a). Using the results presented in Ref. Cvetkovic and Vafek, 2013 and assuming sharply localized Wannier functions of the $xz$, $yz$ and $xy$ orbitals, we can relate, up to a constant, the $d^{\alpha}_{\mu,\sigma}(\textbf{k})$ operators acting in the one iron unit cell to the components of the doublets $\Psi_{X,\sigma}(\textbf{k})$ and $\Psi_{Y,\sigma}(\textbf{k})$ in the two iron unit cell via

[TABLE]

where we absorb the constant prefactors into the Hubbard and Hund terms. As the $xy$ orbital contributes to both X- and Y-point eq.(21) leads to ”Umklapp” terms at the M-point.

We assume that even parity solutions are still the leading pairing instabilities (for odd parity solutions see Appendix B) and use Eqs.(19-21) to project eq.(18) onto the low energy model decoupled into the spin singlet $A_{1g}$ $s$-wave and $B_{2g}$ $d$-wave symmetry states which at the M-point in presence of inter-band SOC couple to the $E_{g}$ even parity spin triplet state. When defining the two doublets $\Psi^{T}_{1\sigma}(\textbf{k})=(c_{yz\sigma}(\textbf{k}),c_{xz\sigma}(\textbf{k}))$ and $\Psi^{T}_{3\sigma}(\textbf{k})=(c_{xy^{X}\sigma}(\textbf{k}),c_{xy^{Y}\sigma}(\textbf{k}))$ the pairing terms read

[TABLE]

The 2x2 block in eq.(27) corresponds to $A_{1g}$ $s$-wave pairing. We write $J^{\prime}_{13}=\alpha J^{\prime}_{11}$ and find two eigenvalues $E_{A_{1g}}=\frac{1}{4}\left(J_{11}^{\prime}2+3U\pm\sqrt{(J^{\prime}_{11})^{2}+\alpha 16J_{11}^{\prime}-2J^{\prime}_{11}U+U^{2}}\right)$ which correspond to ordinary ”plus-plus” ($s^{++}$) $s$-wave and sign-changing ”plus-minus” ($s^{\pm}$) $s$-wave pairing, respectively. While the former channel is purely repulsive without spin-orbit coupling the latter becomes attractive once $J^{\prime}_{11}>(U+U\sqrt{1+8\alpha})/4\alpha$.

The paring term that leads to $B_{2g}$ $d$-wave is $E_{B_{2g}}=\frac{1}{2}(U-J_{11}^{\prime})$ and can be directly read off. It is attractive once $J_{11}^{\prime}>U$ and competes with sign-changing s-wave. Since we assume sharply located Wannier functions which yields eqs.(19-21) and on-site interactions only we find that within our simple mean field approximation the $xy$ orbitals do not contribute to $d$-wave pairing as no ”pair-hopping” term $J^{\prime}_{33}$ mediates between $xy^{X}$ and $xy^{Y}$. This changes once the higher-order diagrams (spin fluctuations) are taken into account.

The $E_{g}$ even parity spin triplet corresponds to pairing between the first (second) component of $\Psi_{X\sigma}(\textbf{k})$ ($\Psi_{Y\sigma}(\textbf{k})$) and the second (first) component of $\Psi_{Y\sigma}(\textbf{k})$ ($\Psi_{X\sigma}(\textbf{k})$) and is thus inter-band and attractive once $E_{E_{g}}=U_{13}^{\prime}-J_{13}<0$. We perform a mean-field decoupling of Eq.(27) into $A_{1g}$ and $B_{2g}$ spin singlet channels with the pairing terms given by

[TABLE]

In terms of $\Psi_{M,\sigma}(\textbf{k})$ a triplet term can be written as $\Psi_{M,\sigma}^{T}(-\textbf{k})\hat{\mathcal{M}}i\sigma^{y}\boldsymbol{\sigma}\Psi_{M,\sigma}(\textbf{k})$, where $\hat{\mathcal{M}}$ and $i\sigma^{y}\boldsymbol{\sigma}$ represent orbital and spin part, respectively. Since $i\sigma^{y}\boldsymbol{\sigma}$ is symmetric an even (odd) parity triplet requires $\hat{\mathcal{M}}$ to be anti-symmetric (symmetric). $i\sigma^{y}\boldsymbol{\sigma}$ can be divided into an in plane $i\sigma^{y}(\sigma^{x},\sigma^{y})$ and out of plane $i\sigma^{y}\sigma^{z}$ component which transform as as the two dimensional $E_{g}$ and one dimensional $A_{2g}$ irreducible representation, respectively. In presence of SOC, orbital and spin degrees of freedom transform together under operations of the space group. We focus on the $E_{g}$ even parity spin triplet that together with the $E_{g}$ in-plane spin component decomposes into a direct sum of one dimensional representations as $E_{g}\otimes E_{g}=A_{1g}\oplus B_{1g}\oplus B_{2g}\oplus A_{2g}$. Using the two anti symmetric components of $E_{g}$ and $\tau_{\pm}=(\tau_{1}\pm i\tau_{2})/2$

[TABLE]

one finds two even parity spin triplets that transform according to $A_{1g}$ and $B_{2g}$ and hence, couple to the singlet channel.

[TABLE]

We refer the reader to Refs.Kang and Fernandes (2016); Eugenio and Vafek (2018) for further details.

In terms of the spinor $\Psi^{T}_{g}(\textbf{k})=(\Psi_{X\uparrow}(\textbf{k}),\Psi_{Y\downarrow}(\textbf{k}),\Psi^{\dagger}_{X\downarrow}(\textbf{k}),-\Psi^{\dagger}_{Y\uparrow}(\textbf{k}))$ the BdG-Hamiltonian reads

[TABLE]

with

[TABLE]

and the pairing terms

[TABLE]

where the gaps in orbital space are given by the equations (LABEL:Eq:s-wave) and (LABEL:Eq:d-wave). We self-consistently compute $\hat{\Delta}_{A_{1g}}$ and $\hat{\Delta}_{B_{2g}}$ as a function of temperature and inter-band SOC for the two cases $E_{A1g}<E_{B2g}$ and $E_{A1g}>E_{B2g}$ and present the phase diagrams in Fig.2(a) and Fig.2(b), respectively. Figure2(c) shows the angular dependence of the superconducting gap projected onto inner (blue) and outer (red) Fermi surface for three values of $\lambda_{\text{SOC}}$ marked by black arrows. The gap projected on the l’th band is given by $\Delta_{l}(\textbf{k})=\sum_{\alpha,\beta}\mu^{\dagger}_{l\alpha}(\textbf{k})\hat{\Delta}_{\alpha\beta}\mu_{\beta l}(\textbf{k})$ where $\hat{\mu}(\textbf{k})$ diagonalizes $\mathcal{H}_{0}(\textbf{k})$.

In the case when solely $E_{A_{1g}}$ ($E_{B_{2g}}$) is attractive, $\Delta^{A}_{t}$ ($\Delta^{B}_{t}$) is induced by SOC and scales with $\lambda_{\text{SOC}}$. In this case $|\Delta^{A(B)}_{t}|<|\Delta^{A(B)}_{1}|,|\Delta^{A(B)}_{3}|$ and we call the state a spin singlet dominated $A_{1g}$ ($B_{2g}$). If $E_{E_{g}}$ is attractive, however, $\Delta_{t}$ can develop independently of $\Delta_{1,3}$ and allows for states where $|\Delta^{A(B)}_{t}|>|\Delta^{A(B)}_{1}|,|\Delta^{A(B)}_{3}|$ which we call spin triplet dominated $A_{1g}$ ($B_{2g}$). Equally whether the system is in the singlet or triplet dominated regime, the form of Eqs. (51) and (56) stays the same.

Depending on the ratio $\alpha=J^{\prime}_{11}/J^{\prime}_{13}$ the leading pairing symmetry for low values of $\lambda_{\text{SOC}}$ is either $B_{2g}$ $d$-wave or $A_{1g}$ $s$-wave with accidental nodes and a sign change between inner and outer Fermi pocket. In Fig.2(a) we show a SOC mediated transition from $B_{2g}$ $d$-wave to $A_{1g}$ s-wave. For small values of $\lambda_{\text{SOC}}$ we find a nodeless $d$-wave state where the nodes are lifted and the gap is opened by a combination of spin singlet and triplet components with the singlet gap being dominant $|\Delta^{B}_{1}|>|\Delta^{B}_{t}|$. With increasing $\lambda_{\text{SOC}}$ the paring symmetry changes from nodeless to nodal d-wave and finally to nodeless s-wave with a dominant inter-band triplet component $|\Delta^{A}_{t}|>|\Delta^{A}_{1}|,|\Delta^{A}_{3}|$. In Fig.2(b) we choose parameters such that the initial state at small $\lambda_{\text{SOC}}$ is $A_{1g}$ $s$-wave and has a dominant singlet gap $|\Delta^{A}_{1}|>|\Delta^{A}_{t}|$. Increasing $\lambda_{\text{SOC}}$ lifts the nodes at the Fermi level and drives the system in an $A_{1g}$ $s$-wave symmetry with a dominant triplet gap $|\Delta^{A}_{t}|>|\Delta^{A}_{1}|,|\Delta^{A}_{3}|$. We find that for large inter-band SOC the triplet dominated $s$-wave state wins over triplet dominated $d$-wave. The presence of an intra-band SOC term $\lambda_{z}$ does not change the phase diagram qualitatively. However, in the region where $\lambda_{\text{SOC}}$ is strong superconductivity in the $d$-wave channel is suppressed if $\lambda_{z}>0$ which further stabilizes the $A_{1g}$-wave solution. Moreover, a large triplet gap, in addition to SOC, lifts accidental nodes leading to a nodeless $s^{\pm}$ pairing symmetry.

Figure2(a) is consistent with the phase diagram of phenomenological model in Ref. Khodas and Chubukov, 2012 with several important differences. In particular, in our model the singlet dominated d-wave state competes with singlet dominated bonding anti-bonding $s^{\pm}$ which arises from pair hopping between $xz(yz)$ and $xy$ orbitals. Furthermore, with increasing spin-orbit coupling strength the $A_{1g}$ state contains dominant interband spin triplet component which was absent in the simplified analysis of Ref. Khodas and Chubukov, 2012 as its pairing channel turned out to be strongly repulsive in the two band model. Since SOC couples even parity inter-band triplet to intra-band even parity spin singlet pairing, an attraction in the former induces a gap in the latter channel and vice versa. Consequently a gap at the FS opens even though the $E_{g}$ triplet is inter-band. A similar paring state with attraction in the triplet channel was recently proposed for highly doped systems with only hole pocketsVafek and Chubukov (2017).

IV Quasiparticle interference

IV.1 Local density of states in presence of impurities

In order to investigate the QPI signatures of the possible pairing states we need to calculate the local density of states (LDOS) for a multi-orbital system in the superconducting state. We further investigate corrections to the LDOS which arise from scattering at single charge impurities to make a statement whether there is a sign-change between inner and outer electron pocket using HAEM’s method. Moreover, we consider scattering from a single magnetic impurity to reveal information about the spin structure in the system.

In order to describe superconductivity and spin-resolved STM we introduce the 16 component Balian-Werthammer spinor

[TABLE]

Within this basis the Hamiltonian in the superconducting state is given by

[TABLE]

where the additional factor of $\frac{1}{2}$ accounts for double counting and

[TABLE]

From Eq.(60) we find the superconducting Green’s function $\hat{G}^{0}_{\textbf{k}}(\omega)=[(\omega+i\delta)\mathbbm{1}_{16\times 16}-H_{BdG}(\textbf{k})]^{-1}$ and the local density of states (LDOS) given by

[TABLE]

where we assumed sharply localized Wannier functions of the $xz$, $yz$ and $xy$ orbitals.

We now introduce a single on-site non-magnetic potential scatterer to the system which is located on either sub-lattice A or B. In the one-iron unit cell it can be described as

[TABLE]

We assume that the major contribution to scattering is of intra-orbital nature and project eq.(67) onto the states in the two iron unit cell. Due to the fact that $d^{A/B}_{xy}$ in the projected model is a linear combination of $c_{xy^{X}}$ and $c_{xy^{Y}}$ the intra-orbital term induces scattering between $xy^{X}$ and $xy^{Y}$ thus also contributes to inter-band scattering. In the two iron unit cell the impurity potential is given by

[TABLE]

with

[TABLE]

Following the HAEMHirschfeld et al. (2015) approach, we compute the antisymmetrized correction to the LDOS due to impurity scattering

[TABLE]

with $\delta\hat{G}_{\textbf{q}}(\omega)=\sum_{\textbf{k}}\hat{G}^{0}_{\textbf{k}}(\omega)\hat{U}\hat{G}^{0}_{\textbf{k}+\textbf{q}}(\omega)$ being the convolution of the bare Green’s functions dressed by a Nambu scattering matrix $\hat{U}=\tau_{3}\otimes\sigma_{0}\otimes\tilde{V}$ in Born approximation, i.e. with $V_{\alpha,\beta}\ll E_{F}$, where $E_{F}$ is the Fermi energy.

The HAEM method states that the momentum integrated and antisymmetrized LDOS, $\delta\rho^{-}(\omega)$, qualitatively differs in case of a sign-changing OP from that of a sign-preserving one yielding a strong enhancement of the integrated response for the sign-changing but not the sign-preserving one. This method has been recently successfully applied to confirm the sign-changing nature of the order parameter in FeSeSprau et al. (2017), where the superconducting gaps are also extremely anisotropic. More recently, the method has also been applied to (Li1-xFex)OHFe1-yZnySeDu et al. (2017) where $\delta\rho^{-}(\omega)$ shows a strong signal and no sign-change between 8 meV and 14 meV suggesting a sign-changing $s^{\pm}$ pairing symmetry. Since our previous analysis shows that most of the conclusions regarding the phase structure of the superconducting gap obtained within Born limit are robust and remain valid also well beyond this limit, we restrict our analysis to weak potential scatterers.

IV.2 Local density of states and phase sensitive correction to QPI

We would like to discuss singlet and triplet dominated $A_{1g}$- and $B_{2g}$-wave with respect to their consistency with QPI experiments in (Li1-xFex)OHFe1-yZnySe. Recently in Ref.(Eugenio and Vafek, 2018) possible pairing states for intercalated FeSe have been discussed. Following Ref.[Eugenio and Vafek, 2018] the pairing has to obey three criteria in order to be consistent with experiments: i) fully gaped LDOS, ii) the quasiparticle energy extrema are at or close to $k_{F}$ of the normal state (”back bending”) and iii) two peak features in the DOS with a peak at about 8 meV and 14 meV, respectively. In the context of monolayer and intercalated FeSe it has been suggested that the second peak in the DOS arises due to pure interband gapEugenio and Vafek (2018) as in the former the gap size seen by ARPESZhang et al. (2016) (13.7 meV) significantly deviates from the peak energy (20.1 meV) seen by STMWang et al. (2012). This, however, does not have to be the case with the electron-intercalated materials where a gap of $13\pm 2$ meV around the Fermi energy was reportedZhao et al. (2016). Three pairing states were proposed to be consistent with experimental data. These are $A_{1g}$ $s$-wave and $B_{2g}$ $d$-wave and $E_{u}\otimes U(1)$ helical $p$-wave. The $A_{1g}$ and $B_{2g}$ state were assumed to have a dominant intra-band singlet gap in order to ensure a full gap and back bending. An inter-band SOC then mixes in inter-band triplet $E_{g}$ pairing. For the odd parity intra-band triplet $E_{u}\otimes U(1)$ $p$-wave state its the odd parity inter-band singlet $A_{2u}$ that is coupled via SOC.

In the following we investigate spin singlet and spin triplet dominated $A_{1g}$- and $B_{2g}$-wave for their consistency with the citeria i)-iii) and the QPI data. For that we use Eq.(60) and calculate the LDOS ($\rho(\omega)$) and the superconducting band dispersion to ensure i)-iii) are fulfilled. In addition, for each pairing state, we calculate the antisymmetrized correction to the LDOS ($\delta\rho^{-}(\omega)$) to check whether the order parameter changes sign between electron pockets. Since the gaps in the intercalated FeSe are found to be isotropicZhao et al. (2016) we further present the angular dependence of the gap projected on the inner and outer electron pocket, respectively. In Appendix (B) we also briefly discuss the QPI data for the odd parity $E_{u}$-wave state and show that at least within Born scattering the results are not compatible with experiment.

We find that $A_{1g}$ and $B_{2g}$ states with a dominant spin triplet gap can show back bending and, in contrast to their singlet dominated states and the odd parity $E_{u}$-state, give an appropriate description of the QPI data found in experiment.

One should also further note that for the HAEM method the crucial role is played by the gaps present on the Fermi surface pockets. In the present case there is also additional interband gap. We show in the Appendix C that its phase structure with respect to the gaps present on the Fermi level cannot be elucidated within phase-sensitive QPI analysis. Furthermore, we show that the sign-changing and sign-preserving gaps on the Fermi surface still determine the characteristic features of $\delta\rho^{-}(\omega)$.

IV.3 A1g-symmetry state

We start by examining spin singlet and triplet dominated $s$-wave pairing state. We use the singlet gaps $\Delta_{1}$ and $\Delta_{3}$ and the triplet gap $\Delta_{t}$ to fit band-structure and LDOS.

IV.3.1 Spin-singlet dominated $A_{1g}$-wave state

To start with singlet dominated $s$-wave we utilize the fitting parameters presented in Ref. Eugenio and Vafek, 2018. As shown in Fig.3(a), $\rho(\omega)$ is indeed fully gaped and has a peak at $8$ meV and $14$ meV. The first peak between $6-10$ meV mostly comes from the intra-band gaps at the electron pockets as can be seen by comparison with the gap projections on the inner, $\Delta_{in}$ and outer, $\Delta_{out}$, electron pocket, Fig.3(c). At $\theta_{F}=\pi/4$ both gaps are equal in magnitude leading to a peak at $\omega\sim 8$ meV in $\rho(\omega)$. The width of the peak is limited by 6 meV and 10 meV which correspond to the minimum of $\Delta_{in}$ and the maximum of $\Delta_{out}$, respectively, so that the occupied states below and above 8 meV are due to anisotropy of the gaps. The second peak at 14 meV is indeed due to inter-band pairing. It vanishes for $\lambda_{\text{SOC}}\rightarrow 0$ and $\Delta_{t}\rightarrow 0$. We find that both gaps are highly anisotropic and most important positive on both Fermi surface pockets, yielding an $s^{++}$-pairing symmetry. In Fig.3(b) $\delta\rho^{-}(\omega)$ exhibits a sign change in the region $6\ \text{meV}<\omega<10\ \text{meV}$ which can be interpreted as a sign preserving order parameter and thus is not compatible with the QPI dataDu et al. (2017). Note that this state does not appear in our phase diagram, shown in Fig.2, where the A1g state is sign-changing. The feature at $\omega\sim 14$ meV in $\delta\rho^{-}(\omega)$ for this state appears due to the small inter-band contribution and does not carry phase information as we explain for a simple model in Appendix C. In Fig.3(d) the superconducting band dispersion in $\Gamma M$-direction ($\theta_{F}=\pi/4$) is shown where the minima of the lower superconducting band are located above the former Fermi level (back-bending).

It turns out that the fitting parameters for spin singlet dominated $s$-wave, proposed in Ref. Eugenio and Vafek, 2018 would give rise to $s^{++}$ superconductivity. A possible $s^{\pm}$ pairing symmetry with a dominant spin singlet gap, on the contrary, would require $\text{sign}\Delta^{A}_{1}\neq\text{sign}\Delta^{A}_{3}$. This state in our anaylsis exhibits highly anisotropic gaps and possible accidental nodes (see for example Fig.2c(1)). This makes it difficult to fit the U-shaped two peaked LDOS in the LiOH-intercalated FeSe data without invoking a spin triplet component, induced by spin-orbit coupling that lifts the nodes.

IV.3.2 Spin-triplet dominated $A_{1g}$-wave state

In Fig.4(a) and Fig.4(d) the LDOS and the superconducting band dispersion are shown. The LDOS shows a full gap and a peak at 8 meV and 13.7 meV, respectively, while the lower superconducting band shows back bending. The band projected gaps $\Delta_{in}$ and $\Delta_{out}$ we plot in Fig.4(c). Possible accidental nodes are lifted by a large $\Delta_{t}$ leading to a nodeless $s^{\pm}$ pairing symmetry with two almost isotropic gaps. In Fig.4(b) we present $\delta\rho^{-}(\omega)$ which in contrast to the singlet dominated case does not change sign between 8 meV and 14 meV. If we compare Fig.4(a) and Fig.4(c) we find that the peak positions of the two peaks roughly agree with the the magnitudes of $\Delta_{in}$ and $\Delta_{out}$ and thus $\delta\rho^{-}(\omega)$ and the behavior of $\delta\rho^{-}(\omega)$ agrees with the experimental one. The first peak is sharper than the second since $\Delta_{in}$ is very isotropic. The second peak is broader as a consequence of an anisotropic $\Delta_{out}$ and the contribution from inter-band pairing. In contrast to the singlet dominated case $\delta\rho^{-}(\omega)$ exhibits well pronounced negative peaks at $\omega_{1}\approx\Delta_{out}$ and $\omega_{2}\approx\Delta_{in}$ without a sign change between them indicating an OP that changes sign between the electron pockets consistent with the experimental QPI data and ARPES data Zhao et al. (2016); Du et al. (2017). Also note that for the spin triplet dominated $A_{1g}$ state the direct gap feature caused by $\Delta_{out}$ is very close to the interband feature in the DOS which agrees with the findings in Ref.Zhao et al., 2016 and Du et al., 2017.

IV.4 $B_{2g}$-symmetry state

We now discuss the singlet and triplet dominated $B_{2g}$ $d$-wave pairing state. We use the singlet gaps $\Delta_{1}$ and $\Delta_{3}$ and the triplet gap $\Delta_{t}$ to fit band structure and LDOS. Note that even though our simple mean field approach leads to $\Delta_{3}=0$ for the $d$-wave case, as it involves decoupling of the repulsive intraorbital Hubbard interaction, a nonzero $\Delta_{3}$ is necessary to achieve qualitative agreement when fitting the experimental data. In theoretical calculations, a non-zero $\Delta_{3}$ appears due to inclusion of the spin fluctuation diagrams in the Cooper-pairing channel yielding momentum dependent interaction, see Ref. Maier et al., 2011.

IV.4.1 Spin-singlet dominated $B_{2g}$-state

For the spin singlet dominated $B_{2g}$-wave state appropriate fitting parameters were found in Ref.(Eugenio and Vafek, 2018). As shown in Fig.5(a) and Fig.5(d) the LDOS is fully gaped and has a peak at $8$ meV and $14$ meV, moreover, the lower superconducting band shows back bending. $\Delta_{in}$ and $\Delta_{out}$ are plotted in Fig.5(c). Both gaps are almost equal in magnitude and nodeless due to smallness of hybridization between the electron pockets, i.e $|\Delta_{in/out}|>\lambda_{\text{SOC}}$. The first peak in the LDOS is due to the intra-band gaps $\Delta_{in}$ and $\Delta_{out}$. We calculate $\delta\rho^{-}(\omega)$ and present the results in Fig.5(b). $\delta\rho^{-}(\omega)$ has positive peaks at $\sim 7$ mev and $\sim 9$ meV with no sign change between them indicating the sign change between $\Delta_{in}$ and $\Delta_{out}$. Yet the third peak at $\sim 14$ meV comes from inter-band pairing and causes two sign changes between second and third peak which is at odds with the experimental dataDu et al. (2017).

IV.4.2 Spin-triplet dominated $B_{2g}$-state

Even though a triplet dominated $d$-wave scenario is not realized in the phase diagram we quickly discuss this case with respect to criteria i)-iii) and QPI data.

As shown in Fig.6(a) and Fig.6(d) the LDOS is fully gaped and has a peak at $8$ meV and $14$ meV, moreover, the lower superconducting band shows back bending. $\Delta_{in}$ and $\Delta_{out}$ are plotted in Fig.6(c). Both gaps are equal in magnitude and nodeless due to smallness of hybridization between the electron pockets. In contrast to the singlet dominated case for the triplet dominated case we use $|\Delta^{B}_{t}|>|\Delta^{B}_{1/3}|$ and $\text{sign}(\Delta^{B}_{t})=\text{sign}(\Delta^{B}_{1})=\text{sign}(\Delta^{B}_{2})$. In Fig.6(b) we plot $\delta\rho^{-}(\omega)$. Two small peaks between 7 meV and 8 meV reflect $\Delta_{in}$ and $\Delta_{out}$ contributions to $\delta\rho^{-}(\omega)$. The peak at 15 meV is due to the large inter-band gap. Overall, the behavior of $\delta\rho^{-}(\omega)$ is nearly consistent with experimental data as $\delta\rho^{-}(\omega)$ does not change sign in the experimentally relevant energy range. However, one has to bear in mind that in contrast to the triplet driven $A_{1g}$ case the phase structure of the gaps on the Fermi surface affect the behavior of $\delta\rho^{-}(\omega)$ only near the lower peak, and the structure near the second (interband) peak is determined by the interband gap and does not necessarily bear information about the signs of the order parameters (see Appendix C.2).

As we clearly see, to obtain the agreement with QPI experimentsDu et al. (2017) the spin-triplet interband component from $E_{g}$ state is necessarily required in both symmetry states and triplet dominated $A_{1g}$ state appears most likely one. This poses an important question how to detect it. To stay within QPI we propose to employ spin-resolved STM as an additional tool to further specify the underlying pairing state and also to distinguish between $B_{2g}$ and $A_{1g}$-states.

IV.5 Spin resolved STM

In the previous section we investigated corrections to the LDOS which arise from scattering on a single non-magnetic charge impurity and found the triplet dominated $A_{1g}$ and $B_{2g}$ pairing states in the presence of SOC to be consistent with the available experiments. As discussed in the introduction, however, the more conventional spin singlet $A_{1g}$ and $B_{2g}$ states have also been proposed. This raises the question whether SOC effects and the spin triplet order parameter can be further verified in experiment. Here we propose spin-resolved QPI as an additional tool to further specify the underlying pairing state. To illustrate the possible capabilities of this technique, we study the corrections to the LDOS in Born approximation from a single magnetic impurity which, within our basis of eq.(57), is given by $\hat{U}=\tau_{3}\otimes\sigma^{z}\otimes\hat{V}$. Here $\tau$ and $\sigma$ matrices act on Nambu and spin-space, respectively, and the matrix $\hat{V}$ describes orbital scattering. The Fourier transform of the spin resolved $\sigma^{i}$-projected correction to the LDOS is given by [Hofmann et al., 2013]

[TABLE]

with the real and the imaginary part being even and odd in q, respectively. Note that a non-vanishing imaginary part requires the scattering potential to break inversion symmetry, which, as we show below can result from interorbital scattering. The scattering matrix $\hat{V}=\hat{V}_{z}+\hat{V}_{x}$ can be written as

[TABLE]

where, in Born approximation, $\hat{V}_{z}$ potentially contributes to the $\sigma^{z}$- and $\hat{V}_{x}$ to the $\sigma^{x(y)}$-polarized LDOS when taking the trace in eq.(76). The labels intra and inter denote intra- and interorbital scattering, respectively. Note that in each block of the matrices $\hat{V}_{z}$ ($\hat{V}_{x}$) it is the off-diagonal (diagonal) elements that break inversion symmetry since under inversion $\Psi_{X}(-\textbf{k})=\sigma^{z}\Psi_{X}(\textbf{k})$ and $\Psi_{Y}(-\textbf{k})=-\sigma^{z}\Psi_{Y}(\textbf{k})$. Those elements contribute to $\text{Im}\delta\rho_{\sigma^{i}}(\textbf{q},\omega)$, while inversion symmetric ones contribute to $\text{Re}\delta\rho_{\sigma^{i}}(\textbf{q},\omega)$.

Results are presented in Fig.7 exemplary for $\omega=8$ meV where we show real and imaginary part of $\delta\rho_{\sigma^{i}}(\textbf{q},\omega)$ as a function of the scattering vector q for a $\sigma^{z}$ impurity with $i=z$ (a-d) and $i=x$ (e-h). Note that a similar behavior of these quantities is found for all frequencies within the range of the two maxima in the LDOS at 8meV and 14meV, respectively.’

In Fig.7.(a,b) $\text{Re}\delta\rho_{\sigma^{z}}$ shows an even in q $C_{4}$-symmetric QPI pattern qualitatively similar for the $A_{1g}$ and $B_{2g}$ pairing state. In contrast to that, see Fig.7.(c,d), $\text{Im}\delta\rho_{\sigma^{z}}$ is odd in q and exhibit a few qualitative differences between $A_{1g}$ and $B_{2g}$ especially along the intensity edges. Note, however, that $\text{Im}\delta\rho_{\sigma^{z}}$ is caused by the inter-orbital terms in $\hat{V}_{z}$ and hence is expected to show a weakened intensity if $V^{\text{inter}}\ll V^{\text{intra}}$.

In Fig.7.(e-h) we present real and imaginary part of the $\sigma^{x}$-polarized QPI pattern $\delta\rho_{\sigma^{x}}(\textbf{q},\omega)$. We would like to stress that this quantity is non-zero only in presence of either a triplet gap $\Delta^{A(B)}_{t}$ or SOC, which then itself induces a triplet gap. The same is true for the $\sigma^{y}$- polarized state which is related to $\delta\rho_{\sigma^{x}}(\textbf{q},\omega)$ by $90^{\circ}$ rotation.

In Fig.7(e-f) we show the real part $\delta\rho_{\sigma^{x}}$ which is $C_{2}$-symmetric and even in q while the imaginary part is odd, see Fig.7(g,h). In contrast to $\delta\rho_{\sigma^{z}}$ the patterns for $A_{1g}$ and $B_{2g}$ representation show qualitative differences in both real and imaginary part. Moreover, if the impurity occupies one Fe site (either on sub lattice A or B) $\hat{V}_{x}$ contains scattering between $xy^{X}$ and $xy^{Y}$ components which breaks inversion symmetry and according to eq.(21) is of intra-orbital nature. Hence, we expect the QPI signal of $\text{Im}\delta\rho_{\sigma^{x}}$ to be of the same order of magnitude as $\text{Re}\delta\rho_{\sigma^{z}}$ and hence may provide a quantity that allows to distinguish triplet driven $A_{1g}$ from $B_{2g}$.

Thus we have shown that the spin-resolved QPI is a useful tool to study the details of the spin-orbit coupling and the triplet order parameter; our results suggest that a combination of spin-resolved QPI with the more conventional probes can be used to determine the symmetry as well as spin structure of pairing in iron-based superconductors.

V Conclusion

In conclusion, we remind the reader that the highly electron-doped Fe-chalcogenides, in particular the (Li1-xFex)OHFeSe system, have been discussed intensively in the context of the standard model of spin-fluctuation induced spin singlet pairing in the limit of weak spin-orbit coupling. In this case, the nodeless $d$-wave ($B_{2g}$) and bonding-antibonding $s_{\pm}$ ($A_{1g}$) spin singlet states on the electron pockets have been considered most favorable, with the additional possibility of incipient $s_{\pm}$ pairing via coupling to the incipient hole pocket at $\Gamma$. One of these states may eventually prove to be the correct pairing state for these materials, if the effective attraction in these channels generated by spin fluctuations driven by interpocket repulsion dominate the pair vertexMaier et al. (2011).

Here, we have instead studied the new possibility raised previouslyCvetkovic and Vafek (2013), namely that interorbital triplet components generated by intrinsic attractions possible only in the presence of spin orbit coupling are responsible for the pairing in these systems. We explored the various possible Cooper-pairing symmetries using a mean-field decomposition of the Hubbard-Kanamori Hamiltonian including $A_{1g}$ and $B_{2g}$ states. For nonzero spin-orbit coupling, the superconducting order parameter is a combination of spin singlet and spin triplet gaps in each state. Treating attraction in singlet and triplet channels on equal footing and solving the selconsitency equations, we found for weak spin-orbit coupling a dominant spin singlet and small spin triplet gap yielding a state essentially equivalent to those identified in the usual spin fluctuation approach. For stronger spin-orbit coupling, however, the superconducting order parameter is a combination of spin singlet and dominant spin triplet gaps in each state. Focusing on the (Li1-xFex)OHFeSe system, we identified the even parity $A_{1g}$- and $B_{2g}$- pairing states with a dominant spin triplet component to be consistent with available experiments, including current quasiparticle interference data, whereby according to our phase-diagram the $A_{1g}$ state is slightly favored. The spin-singlet dominated $A_{1g}$ and $B_{2g}$-states in this scenario without strong spin fluctuations are not consistent with at least one of the existing experiments.

In summary, to obtain a full moderately anisotropic gap on the Fermi pockets and its sign-changing character in agreement with experimental results on (Li1-xFex)OHFeSe, we require either the traditional intraband $A_{1g}$ or $B_{2g}$ states generated by spin fluctuations, or the new triplet interband pair states in the same symmetry channels generated by intrinsic attraction in multiorbital correlated models. It is clearly of interest to identify experimental tests to distinguish between these possibilities. To this end we proposed using spin-polarized QPI to identify the possible triplet components present in the more exotic alternative states, and presented results for each of the triplet dominated states. A clear identification of triplet interband pairing using these results would be an important step forward in understanding the unusual superconductivity in the Fe chalcogenides.

Acknowledgements.

We thank P. M. Eugenio, O. Vafek, A. V. Chubukov and P. Coleman for useful discussions. J.B. and I.E. were supported by the joint DFG-ANR Project (ER 463/8-1) and DAAD PPP USA N57316180. P.J.H. was supported by the U.S. Department of Energy under Grant No. DE-FG02-05ER46236. P.A.V. acknowledges the support by the Rutgers University Center for Materials Theory Postdoctoral fellowship.

Appendix A Mean-field order parameters

Using mean-field we decouple eq.(27) into $A_{1g}$, $B_{2g}$ and $E_{g}$ pairing channels. The pairing terms for $A_{1g}$- and $B_{2g}$-wave states read

[TABLE]

[TABLE]

respectively.

Appendix B $E_{u}$-wave state

The odd parity $E_{u}\otimes U(1)$ $p$-wave is a novel state which was proposed for intercalated FeSe in Ref. Eugenio and Vafek, 2018. It is a time reversal symmetric topological superconductor which is nodeless for non zero intra-band SOC $\lambda_{z}$. Consequently the order parameter is sensitive to the changes in $\lambda_{z}$ and the values for $\lambda_{z}$, $p_{z_{1}}$ and $p_{z_{2}}$ need to be adjusted to obey criteria (i)-(iii). These were found in Ref. Eugenio and Vafek, 2018 and are $\lambda_{z}=-31\ \text{meV}$, $p_{z_{1}}=4179.88\ \text{meV}\mathring{A}^{3}$ and $p_{z_{2}}=2.5p_{z_{1}}$. This state does not appear in our mean-field analysis. However, for completeness we also calculate the QPI and spin-resolved QPI data for the $E_{u}\otimes U(1)$ helical p-wave symmetry where the pairing term within the basis of eq.(57) reads

[TABLE]

and

[TABLE]

In Fig.8(a) and Fig.8(d) we present the DOS and the energy band dispersion for the $p$-wave state. The DOS show two peaks tuned to lie at 8 meV and 14.3 meV, respectively and the lower energy band shows the back bending. Fig. 8(c) shows the absolute values of $\Delta_{in}$ and $\Delta_{out}$ projected on the Fermi surface pockets. In Fig.8(b) $\delta\rho^{-}(\omega)$ is shown. Here, one finds that the main feature of $E_{u}\otimes U(1)$ in $\delta\rho^{-}(\omega)$ is that there are two sign changes between the peak energies of the total DOS (i.e. two negative peaks and positive values in between them), which does not agree with experimentDu et al. (2017). Therefore this state does not seem consistent with the QPI data at least within Born scattering limitDu et al. (2017). For consistency the spin-resolved QPI for the odd parity $E_{u}\otimes U(1)$ p-state is presented in Fig.9.

Appendix C HAEM’S method for strong inter-band pairing

The SOC induced coupling between the spin singlet and triplet pairing channels translates into a coupling between intra- and inter-band order parameters in band space. It has been arguedEugenio and Vafek (2018) that the second peak, seen in STMDu et al. (2017) can be at least partially due to an interband gap. In order to investigate how the predictions of HAEM’s theory are affected by inter-band pairing we consider a simple model of a two band superconductor with superconductivity driven by intra- and inter-band pairing.

[TABLE]

where $\Psi^{\dagger}(\textbf{k})=\left(c^{\dagger}_{1}(\textbf{k})\enspace c^{\dagger}_{2}(\textbf{k})\enspace c_{1}(-\textbf{k})\enspace c_{2}(-\textbf{k})\right)$. The band dispersions are assumed to be simple parabolic ones $\xi_{1}(\textbf{k})=\frac{k^{2}}{2m}-\mu_{1}$ and $\xi_{2}(\textbf{k})=\frac{k^{2}}{2m}-\mu_{2}$. Hence, $\xi_{2}=\xi_{1}-2B$ with $B=\frac{1}{2}(\mu_{2}-\mu_{1})>0$. We only distinguish between intra- and inter-band pairing and neglect for a moment their spin symmetry (i.e. consider them to be spin singlet) to get a better understanding on how the QPI data is affected by the inter-band pairing. We linearize both bands, and calculate $\delta\rho^{-}(\omega)$ and investigate how it depends on $\Delta_{12}$, the relative phase between $\Delta_{12}$ and $\Delta_{21}$ and the band offset (direct gap).

C.1 QPI for Inter-band

We start by examining $\delta\rho^{-}(\omega)$ for zero intra-band gaps $\Delta_{11}=\Delta_{22}=0$ and assume $|\Delta_{12}|=|\Delta_{21}|$. From the Green’s function $G(\textbf{k},\omega)^{-1}=\left[\omega+i\boldsymbol{\delta}-H_{BdG}\right]$ we calculate the momentum integrated Green’s function as

[TABLE]

where $\boldsymbol{\delta}=\delta\text{sign}(\omega)$. The local LDOS we obtain from $\rho(\omega)=-\frac{\text{sign}(\omega)}{\pi}\sum_{i=1}^{2}\sum_{\textbf{k}}\text{Im}\hat{G}(\textbf{k},\omega)_{ii}$. The result is

[TABLE]

which now explicitly includes the interband pair gaps $\Delta_{12}$ and $\Delta_{21}$. Eq. (114) has four poles.

[TABLE]

For $\omega>0$ one needs to consider the cases: i) $B<|\Delta_{12}|$ and ii) $B>|\Delta_{12}|$ where in the former $\omega_{1}$ and $\omega_{3}$ and in the latter $\omega_{2}$ and $\omega_{4}$ correspond to the peak energies in $\delta\rho^{-}(\omega)$.

i) $B<|\Delta|$

In Fig.10, we plot $\delta\rho^{-}(\omega)$ in arbitrary energy units as a function of $\omega$ and $B<|\Delta|$. For $\Delta_{12}=\Delta_{21}$ the behavior of $\delta\rho^{-}$ is what we define here as “odd”, meaning that between the two intraband gap energies $\delta\rho^{-}$ changes sign, so the QPI-pattern is $s^{++}$-like (solid blue curve). In case of $\Delta_{12}=-\Delta_{21}$ $\delta\rho^{-}$ is “even”, and one obtains as ”$s^{+-}$-like” pattern (solid orange curve). Hence, if the inter-band gap is larger than the band offset HAEM is sensitive to the relative phase between $\Delta_{12}$ and $\Delta_{21}$.

ii) $B>|\Delta|$

If $B>|\Delta_{12}|$ the situation is a different one as can be seen by the dotted curves in Fig.10. In both cases $\Delta_{12}=\pm\Delta_{21}$ and $\delta\rho^{-}$ is “odd”, mimicking an $s^{++}$-pattern. Consequently, HAEM’s method is not sensitive to the relative phase between $\Delta_{12}$ and $\Delta_{21}$. Therefore, in a system with a large dominant inter-band gap the QPI signal depends not only on the relative phase between $\Delta_{12}$ and $\Delta_{21}$ but also on the ratio between $|\Delta|$ and $B$, hence on the band structure.

C.2 Inter+Intra-Band

The influence of the inter-band gap on the QPI does affect the HAEM results on sign-changing and sign preserving intra-band gaps. To show this we numerically present the $\delta\rho^{-}(\omega)$ for $\Delta_{11}\neq\Delta_{22}$ and $\text{sign}(\Delta_{11})=-\text{sign}(\Delta_{22})$ leading to an initial $s^{\pm}$ pattern. Then we increase $|\Delta_{12}|$ and show that at a certain magnitude the phase information is lost. The results are plotted in Fig.11 where in the left and right panel we have $\Delta_{12}=-\Delta_{21}$ and $\Delta_{12}=\Delta_{21}$, respectively.

Our results for the LDOS are presented in Fig.11(a) and Fig. 11(b). If $\Delta_{12}=0$, see orange curve, only intra-band gaps $\Delta_{11}=2$ and $\Delta_{22}=-4$ are present which cause a two peak feature in the LDOS which we call peak 1 and 2. If we switch the inter-band paring $\Delta_{12}=1$ on, see black curves, peak 3 and 4 appear as a consequence of inter-band pairing, whereas peak energies 1 and 2 are only slightly affected. If $\Delta_{12}$ is increased further, see blue curves, peak 2 and 3 merge.

In the absence of inter-band pairing the QPI pattern in Fig.11(c)-(d) is of $s^{\pm}$ type, see orange curve. For moderate values $\Delta_{12}=1$ one can nicely distinguish between intra- and inter-band contributions Peak 1 and 2 show an $s^{\pm}$ pattern due to the sign change between $\Delta_{1}$ and $\Delta_{2}$, whereas between peak 3 and 4 the pattern is an $s^{++}$ since the inter-band gap is smaller that the band offset B, see Fig.10. For moderate values of $\Delta_{12}$, see blue curves, peak energies 1 and 2 start to deviate from $\Delta_{1}$ and $\Delta_{1}$ as marked by the vertical dashed lines. As soon as peak 2 and 3 merge the $s^{\pm}$ pattern between peak 1 and peak 2 becomes more and more difficult to resolve. The blue curve corresponds to the situation where the inter-band gap $\Delta_{12}$ is larger than the intra-band gap but smaller that the band offset. The latter condition causes an $s^{++}$ pattern between peak 3 and 4 independent of the relative phase between $\Delta_{12}$ and $\Delta_{21}$, sec.C.1. This is accompanied by $\Delta_{11}$ and $\Delta_{22}$ having a $s^{\pm}$ symmetry. The combination of both patterns leads to a sign-change of $\delta\rho^{-}(\omega)$ in the region $2<\omega<4$ which might be misinterpreted as a an $s^{++}$ pattern between $\Delta_{11}$ and $\Delta_{22}$ and hence carries no concrete phase information.