Odd-frequency superconductivity in dilute magnetic superconductors

We show that dilute magnetic impurities in a conventional superconductor give origin to an odd-frequency component of superconductivity, manifesting itself in Yu-Shiba-Rusinov bands forming within the bulk superconducting gap. Our results are obtained in a general model solved within the Dynamical Mean Field Theory. By exploiting a disorder analysis and the limit to a single impurity, we are able to provide general expressions that can be used to extract explicitly the odd-frequency superconducting function from scanning tunneling measurements.


I. INTRODUCTION
In superconductors, Fermi statistics imposes that the superconducting pairing function is antisymmetric under exchange of the two electrons forming the Cooper pairs.The pairing function must therefore change sign under the exchange of the quantum numbers labeling the two electrons like position, time, orbital index, spin, etc.In the most conventional single band spin-singlet superconductivity 1 , it is the spin component of the pairing function which is antisymmetric, while the space part is symmetric (e.g.s-wave).In spin-triplet superconductivity, often advocated in the so-called ferromagnetic superconductors 2 , such as UGe 2 , URhGe, and UCoGe, the spin component is symmetric while the space component is antisymmetric (e.g.p-wave).Sign change may occur also in other degrees of freedom.More than 40 years ago, Berezinskii proposed that the antisymmetric contribution to the pairing function may derive from the exchange of the electron time coordinates 3 .He proposed that such a situation occurs in 3 He, where the space (swave) and spin (triplet) components would be symmetric but the time component antisymmetric.In this case, in the space reciprocal to time, the pairing is an odd function of frequency (odd-ω) 4,5 .Since then, odd-ω pairing was predicted to be quite a general phenomenon in superconducting systems 6,7 , including for example disordered superconductors 8 and heavy-fermion superconductors 9 .
In more recent years it has been realized that an oddω pairing component can arise when superconductivity is induced in ferromagnetic systems by proximity with a conventional superconductor 10 .In this case the breaking of time-reversal symmetry induced by an effective magnetic field can change the spin-component of the pairing function from being anti-symmetric to symmetric, favoring the appearance of a time-antisymmetric component.Such systems offer the advantage to be artificially built and controlled, opening a path towards applications in the field of spintronic devices 11,12 .More generally, it has been shown that an odd-ω component can indeed arise whenever symmetry breaking occurs, like e.g. the spatial symmetry breaking in non-magnetic junctions 13,14 .
These phenomena have gained more and more interest with the advent of topological materials, where the competition between superconductivity and magnetic orders is often a key ingredient.For instance, proximity effect on a dense chain or wire of magnetic atoms deposited on top of a superconductor gives rise to unconventional superconductivity, marked by the appearance of Majorana edge states localized at the extremities of the chain [15][16][17][18][19] .These emergent degrees of freedom promise to be the fundamental building blocks in the development of quantum computers 20,21 .Another interesting system is realized by magnetic islands on the surface of conventional superconductors [22][23][24] .Here some experimental signatures of topological superconductivity and chiral Majorana edge channels have been reported 25,26 .Understanding the role played by odd-ω superconductivity in such Majorana systems and finding what are its experimental signatures, are fundamental open questions.
There has been then a remarkable effort from a large part of the condensed matter community to experimentally reveal odd-ω pairing.Proposals include, e.g., measurements of thermopower in superconductor-quantumdot-ferromagnet hybrid systems 27 and Josephson effect in superconductor-ferromagnet junctions 4 .Despite this, the detection of odd-ω pairing has remained a theoretical chimera and only very recently its experimental realization has apparently been confirmed.One recent study has reported odd-ω superconductivity at the interface of a topological insulator with a conventional superconductor 28 .In another recent development, which stands upon earlier theoretical work 29,30 , the presence of an odd-ω component in scanning tunneling spectroscopy (STS) has finally been reported in a system of a magnetic impurity in contact with a conventional superconductor 31 .
On the footsteps of that latter study, the goal of the present work is to go beyond the single impurity system, and show that odd-ω pairing can be induced in a conventional superconductor by the collective effect of a finite concentration of magnetic impurities.For this purpose, we shall consider a general model of magnetic impurity sites embedded within a conventional superconducting lattice (as portrayed in Fig. 1): a dilute magnetic super-

conductor (DMS).
The article is organized as follows.In Sec.II we introduce our model and the Dynamical Mean Field Theory (DMFT 32 ) method, that allows us to solve it in a well controlled infinite dimensional limit.Moreover, DMFT is a mean field theory based on Green's functions, thus has the advantage of providing the local spectral functions, which may be directly observed in STS experiments.Our results are presented in Sec.III.In the first part, we show the appearance of impurity bands inside the superconducting gap, the Yu-Shiba-Rusinov bands, which possess an odd-ω component.Differently from the single impurity case 30,31 , odd-ω superconducting pairing is present for the whole system, for both magnetic and non-magnetic sites, making it eventually exploitable in transport and device making.We develop an impurityconcentration scaling analysis that allows us to derive an explicit expression relating the odd-ω superconducting function to the STS local density of states.This enables us to explicitly extract the odd-ω superconducting pairing function, which we compare with the exact DMFT solution.In the second part, we provide a detailed analytic proof of the relations previously derived by considering the diluted disorder limit of the DMFT solution.Finally, Sec.IV provides a summary of our results in order to motivate future experimental investigations in dilute magnetic superconductors.

II. MODEL AND METHOD
The DMS model portrayed in Fig. 1 is described by the following Hamiltonian Here the operator c † i,α creates an electron with spin projection α at site i and t ij is the hopping amplitude between neighboring sites ij .At random sites l we place classical magnetic moments S θ l (with components θ = x, y, z), which couple to the electrons via an exchange parameter J (here σ θ αβ are the Pauli matrices with the spin indices α, β).The ∆ i term describes superconducting pairing, and we consider ∆ i = ∆ on non-magnetic sites while ∆ i = 0 on magnetic sites.The particle density is fixed by the chemical potential µ, and we fix µ i = µ on non-magnetic sites and µ i = µ + δ µ on magnetic sites, where δ µ is an energy offset that describes potential disorder.Magnetic sites are uniformly distributed, randomly occurring at a given site with probability x (0 ≤ x ≤ 1).It was shown in Ref. 31 (Supplemental Material) that spin-orbit coupling (SOC) does not affect the local Green's functions in the case of a single magnetic impurity in the superconductor.We then expect the effect of SOC on the local Green's functions to be negligible in the case of dilute magnetic impurities, x 1, and do not consider it in the Hamiltonian.
As mentioned in the introduction, this many-body Hamiltonian can be solved in a well-controlled fashion in the infinite dimensional limit, where the J coupling can be treated beyond perturbation theory by means of DMFT 32 .Within this method the interaction reduces to a purely local term, and the full lattice problem is mapped on a quantum impurity model coupled to an effective bath of non-interacting fermions.In this case, we need to consider only the on-site one-particle propagator, which in the superconducting state can be conveniently expressed as a 2×2 Nambu matrix related to the Nambu Here G α (ω) and F (ω) are the normal and anomalous components of the Fourier transform on the Matsubara axis of the Nambu Green's function Ĝ(τ ) = − T ψ(τ )ψ † (0) , where T is the time ordering operator for imaginary time τ .The goal in the DMFT approach is to determine the right hybridization function, or bath, of the quantum impurity, which satisfies the DMFT self-consistency condition, that is specified by the original lattice model of Eq. (1).In our problem, however, there are two inequivalent lattice sites, the magnetic and the nonmagnetic ones.Thus, we must consider two separate quantum impurity problems, one magnetic and the other non-magnetic, which become coupled through the selfconsistency condition.
We introduce some simplifying assumptions, which should not qualitatively change the nature of the model behavior.Firstly, the classical magnetic moments S l are assumed frozen, acting as magnetic spin disorder.With this assumption, the DMFT method turns out to be equivalent to the treatment of magnetic disorder within the so-called coherent potential approximation (CPA) [33][34][35] .Secondly, we assume that the magnetic impurities are fully polarized and keep only the S z component in the Hamiltonian of Eq. (1).Finally, without loss of physical generality, in the infinite dimensional limit it is convenient to adopt a Bethe lattice with hopping t ij = t/ √ z, where z → ∞ is the number of first neighbors of each site, whose density of states is a simple semicircle 32 .This greatly simplifies the DMFT equations, which can be written as the two coupled Green's functions equations at non-magnetic (nm) and magnetic (m) sites The model parameters are physically reasonable, as the conventional superconductor has a relatively small gap and the non-interacting density of states at the chemical potential is featureless.

III. RESULTS
A single magnetic impurity interacting with the electrons in a superconductor gives origin to electron bound states at the impurity site known as Yu-Shiba-Rusinov (YSR) states.This YSR state appears as a sharp resonance within the superconducting gap in the spectral density of states, which can be experimentally revealed, e.g., by STS.Such a phenomenon is reminiscent of the bound states in dilute magnetic semiconductors 36 , where the role of the superconducting gap is played by the semiconductor gap.In analogy with this latter case, when a finite concentration x of magnetic impurities is embedded within a bath of electrons, the impurity electrons can communicate via the bath, giving origin to Shiba bands, which also appear within the superconducting gap.
This is indeed what we find in our calculation, as shown in Fig. 2(a), where we display the electronic density of states N ↑/↓ (ω) = −Im[G m↑/↓ (ω)]/π at magnetic sites in the case with x = 0.01.An odd-ω component also appears in non-magnetic sites by an inverse proximity effect, as the odd-ω superconductivity spreads throughout the system.In those sites, however, the amplitude is much smaller (see Appendix A for an example).Therefore in the following we shall focus on the magnetic sites, where it would be easier to detect the odd-ω contribution.Two YSR impurity bands are clearly visible within the superconducting gap, as displayed in the inset of Fig. 2(a).With this choice of parameters, the lower (ω < 0) and upper (ω > 0) YSR bands are well separated and have spin down (dashed line) and spin up (continuous line) character.We shall discuss later the more complicated case where the lower and upper YSR bands overlap close to ω = 0.A possible experimental realization of this has been discussed in conjunction with the presence of Majorana fermions, which appear as a resonance in the density of states at ω = 0 15,18 .
In contrast to dilute magnetic semiconductors, our system is in a superconducting state.We expect therefore superconductivity to be induced at the impurity site by proximity effect.Because of the time reversal symmetry breaking due to the magnetic field at the impurity, an odd-ω pairing component is expected to appear 4 , similarly to the single impurity case 30,31 .As the impurity electrons form the YSR bands, these should then dis-play an odd-ω superconducting component in the superconducting function.In fact, such component can be seen in Fig. 2, panels (b) and (c), where we display the real and imaginary parts of the full superconducting function, Re[F m (ω)] and Im[F m (ω)], respectively.In a non-magnetic spin-singlet BCS-like superconductor Re[F m (ω)] is symmetric in ω, whereas here it is not, in close analogy to what it is found in the single-impurity case 30,31  While the density of states (Fig. 2(a)) can be directly obtained from STS measurements, it is in general difficult to extract the superconducting function and, in particular, the odd-ω part (Fig. 2(b)-(c)).Recently it was shown that this can be achieved in the case of a single isolated impurity, where the odd-ω Cooper pairs are localized at the impurity site 31 .We shall now show, that this idea can be extended to the present case of delocalized YSR bands with odd-ω superconductivity.We shall provide a explicit protocol to extract the superconducting function from STS experimental data.To this purpose, we shall use the magnetic disorder concentration, x, as a tuning scaling parameter and derive explicit expressions relating the STS density of states to the odd-ω superconducting funtion.

A. YSR-band scaling with impurity concentration
We first notice that the shape of the YSR bands (inset of Fig. 2(a)) reflects the semicircle shape of the Bethe lattice density of states D(ε).The presence of a quasiparticle band, with a renormalized mass and a density of states that reflects the non-interacting one of the lattice, is a well-known feature of the DMFT solution for strongly interacting metallic states 37 .Here, the strongly correlated states are those of the magnetic impurity network, whose states are subject to the local interaction J that reduces their effective hopping.Thus, we may expect that this quasiparticle heavy-band may carry the information of the impurity concentration x.Thus, we may attempt a scaling of the YSR bands as a function of the disorder-site density x.In Fig. 3(a) we plot the upper YSR band for various values of x scaled as where E 0 , b, a ↑ and a ↓ depend on the model parameters t, µ, δ µ , ∆ and J, but do not change significantly with x.Note that ω = ±E 0 mark the midpoints of the YSR bands.Also notice that a ↑/↓ have a priori different values for each of the two Shiba bands.Magnetic disorder acts to rescale the width and the height of the YSR bands according to a √ x dependence.The collapse of all the curves on the same line proves the validity of such a scaling.
A similar scaling applies to the superconducting function as well.In Fig. 3(b) we show that the odd-ω component Im[F o m (ω)] has the same shape of D(ω) and the same scaling with x, These results bring us to establish a useful relation between the odd-ω superconducting function and the density of states This relation is similar to the one obtained for the single impurity case 31 , as we show explitely in the following section III B. However its range of validity is now extended to dirty superconductors with dilute magnetic impurities.Equations ( 4)-( 6) remain valid even when the Shiba bands overlap (see section III B).STS measurements of the spin-polarized density of states could allow us, via Eq.( 6), to extract the odd-ω component of the superconducting function, provided the coefficient a F is determined.This cannot be obtained directly from a spectroscopic measurement.We show however in the following section that a relation between a ↑/↓ and a F can be explicitly derived by going to the single impurity limit, as these coefficients depend only weakly on the impurity concentration x.If the normal-state density of states at the Fermi level can be assumed constant within the superconducting gap (i.e.particle-hole symmetry can be assumed at low energy), which is the case of most of standard superconductors, we shall show that With this relation we find the following expression that relates the odd-frequency superconducting function directly to the spin-polarized density of states N ↑/↓ (ω), which are a priori directly measurable functions: Notice that in most cases, when the lower and upper YSR bands do not overlap significantly at ω = 0, the YSR bands are fully spin polarized (as in Fig. 2(a)).We can in this case identify N ↑ (ω) = N (ω > 0) and N ↓ (ω) = N (ω < 0), where N (ω) = N ↑ (ω) + N ↓ (ω) is the total density of states which is more easily accessible with standard STS.We can then simply write Furthermore, even when the lower and upper YSR bands do overlap at ω = 0, it would be still possible to attempt a fit of the experimental density of states which could separate the N ↑/↓ (ω) components.Then, Eq. ( 8) can still be used as a good approximation if the YSR bands overlap over a short interval inside the superconducting gap.Finally, when even this last case is not applicable, it would be still possible to measure N ↑/↓ (ω) components by employing spin-dependent STS 38 .We display in Fig. 4 the Im[F o m (ω)] extracted using Eq. ( 8) for different values of the model parameters, including a case where the YSR bands overlap close to ω = 0.The curves (dashed lines) are compared with the Im[F o m (ω)] obtained directly from the solution of the DMFT equations (3) (continuous lines).The good agreement is a proof of principle of our analysis and validates our results.We finally remark that once Im[F o m (ω)] is extracted using the method presented here, Re[F o m (ω)] can be obtained via the Kramers-Kronig relations.

B. Derivation of relations
We shall now prove that the relation established from the disorder concentration scaling, Eq. ( 6), can be justified analytically and we derive Eq. ( 7) explicitly by going to the well controlled single impurity limit.
In the small concentration regime x 1, magnetic impurities affect only mildly the bulk superconducting Green's function Ĝs .We can then assume Ĝnm Ĝs .
From this relation we derive the following equations (see Appendix B), where G s (ω), F s (ω) are elements of the matrix Ĝs (ω).These expressions illustrate well the DMFT results, as we show in Fig. 6(b) of Appendix B.
Here we follow Ref.31, generalizing it to the case of many magnetic impurities, to obtain some useful expressions for low energies.For |ω| < |∆|, inside the superconducting gap, we can approximate the denominator of Eq. ( 10) by F s (ω) − ∆det Ĝs (ω) Multiplying this equality by (G s (ω) − G * s (−ω))/[πF s (ω)] and adding the result to Eq. ( 10), we find where Here we assumed the same density of states for both spin species in the clean superconductor.These functions can be simplified for ω inside the gap, where the imaginary components of G s (ω) and F s (ω) vanish.For a given lattice density of states D( ), these coefficients can be expressed in terms of model parameters (see Appendix C), by making a power series expansion of D( ) around = µ and keeping the first order terms.We find Here W is a cutoff proportional to the bandwidth and D ( ) is the derivative of the density of states.The function Im[F m (ω)] vanishes outside its resonance around ω = E 0 (see e.g.Fig. 2(c)) and hence the function Im[F m (−sgn(E 0 )|ω|)] that appears in Eq. ( 12) vanishes for every ω where the YSR bands do not overlap.If x is small enough, these bands are expected to overlap only for |ω| |∆|, where C r (ω) is negligible (see Eq. ( 14)).Therefore the last term of Eq. ( 12) will be neglected.
The function Im[F e m (ω)] is antisymmetric while , C e (ω) and C o (ω) are symmetric with respect to ω. Extracting the symmetric (s) and antisymmetric (a) components of Eq. ( 12) leads to We now notice that the functions ] vanish quickly as one moves away from ω = E 0 , while C e/o (ω) vary slowly near ω = E 0 .Therefore, we can safely replace C e/o (ω) with C e/o (E 0 ): Equation ( 16) is then a good approximation even when the YSR bands overlap with one another, if the overlap takes place over a short interval inside the gap.We now derive the relation between the odd-ω superconducting function and the density of states, Eq. ( 6), that we inferred from the disorder concentration scaling analysis of the YSR bands.We first notice that integrating Eq. ( 4) one finds Plugging these relations into Eq.( 16) one recovers Eq. ( 6) provided that Note that C e (E 0 ) is given by Eq. ( 14) and in the particular case where |E 0 | |∆| one has C e (E 0 ) sgn(∆)2JD(µ).Via the relation (6) we can finally extract the odd-ω superconducting function from a disorder concentration scaling analysis of the YSR bands, if we can determine the magnetic coupling J.As the latter may be difficult to extract from experiments, we can further simplify Eq. ( 17) by going to the single impurity limit.
As we have mentioned in the disorder concentration scaling analysis, a ↑/↓ and a F do not depend on x.They thus remain the same in the single impurity limit, x → 0, and Eq. ( 4) becomes a delta-like Lorentzian, where we introduced a finite inverse lifetime parameter η for the Shiba state, replacing ω → ω + iη.Similarly, Eq. ( 6) can be used together with Eq. ( 18) to find a Lorentzian expression for Im[F o m (ω)]: On the other hand, we can take the DMFT Eqs.(3) in the x → 0 limit (see Appendix D for details).Once again we replace Ĝnm Ĝs , and obtain the Lorentzian form of N ↑/↓ and Im[F o m (ω)] as a function of the model parameters D(ω), J, ∆, µ, δ µ (see Eq. (1)), provided that where we define 1, i.e. we assume the density of states is almost constant at low energies.Then we finally obtain which proves the relation (7).

IV. CONCLUSION
We have shown that odd-ω superconductivity occurs in a dilute magnetic superconductor model, where a finite concentration of magnetic impurities are randomly embedded in a conventional superconductor.We have solved this model by treating the coupling between magnetic impurities and superconducting electrons in Dynamical Mean Field Theory.This method is also appropriate to deal with the disorder and, within our approximations, is equivalent to the coherent potential approximation.This technique allows us to have direct access to the local spectral functions that can be experimentally extracted by scanning tunneling experiments.Our results show the formation of Yu-Shiba-Rusinov bands inside the bulk superconducting gap displaying a relevant odd-frequency component.We have analyzed the YSR bands by means of a scaling analysis and derived an expression (Eq.( 8)) that allows us to explicitly extract the odd-frequency superconducting function from spectroscopic quantities that are directly accesible in scanning tunneling measurements.To get this explicit formula we have exploited the findings of our disorder-concentration scaling by going to the single impurity limit, where exact expressions could be derived within physically reasonable approximations.We benchmarked our approximate expressions and found a very good correspondence with the Dynamical Mean Field Theory results, providing further proof to the validity of our analysis.
Our results should motivate future experimental investigations in dilute magnetic superconductors to search for odd-frequency superconductivity.In particular, the formation of YSR bands possessing odd-frequency character, which delocalize into the bulk superconductor, should raise further interesting questions about the impact of the odd-frequency pairing not only on the spectral properties but also on other thermodynamic and transport properties, with an eye for future spintronic device applications.3), as presented in Fig. 2. In (a), dashed red line corresponds to Eq. (B5), showing that it is exact for the whole spectrum.In (b), dashed red line corresponds to the approximate Eq. ( 10), showing that the latter is a good approximation for the former for a broad range of energies.
Eq. (B1) we obtain the relations From Eq. (B2) one finds an expression for det Ĝm (ω) and using it in Eqs.(B3)-(B4) one finds  10) is a good approximation to the DMFT solution for a broad range of energies.However this approximation might be inappropriate for high energies, |ω| ∼ W , where Ĝnm (ω) may be more strongly affected by the impurities, and the approximation Ĝnm (ω) Ĝs (ω) breaks down.
The coefficients in these equations can be determined as follows.

FIG. 1 .
FIG.1.Dilute Magnetic Superconductors: magnetic impurities sites are embedded in a superconducting lattice.The impurity-site magnetic moment S (represented by the blue arrow) interacts with electrons (in orange) via a magnetic coupling J, as described by the Hamiltonian in Eq. (1).
) Here Ĝav = x Ĝm + (1 − x) Ĝnm is the Green's function averaged over m and nm sites and τ x,z are Pauli matrices in the Nambu spinor indices.These equations are solved numerically.For definiteness, in the following figures we fix the Hamiltonian parameters t = 1, µ = −0.05,δ µ = −0.5, ∆ = −0.1 and J = −0.65.The results that we shall describe next are rather generic, namely, they do not depend on any particular fine tuning of parameters.
. It is convenient to decompose the total F m (ω) in a standard even-ω F e m and an odd-ω F o m components, defined by F e/o m (ω) = Fm(ω)±Fm(−ω) * 2 , which are also displayed in Figs.2(b)-(c).The function F o m (ω) describes triplet, s-wave, odd-ω pairing.

FIG. 3 .
FIG. 3. Plots of Shiba bands centered at E0 > 0 using the same parameters of Fig. 2, but with η = 10 −4 and making a scaling with x.The horizontal axes are multiplied by x −1/2 and the vertical ones by x 1/2 .(a) shows plots of N ↑ using x = 0.0001, 0.001, 0.005, 0.01, 0.015, 0.02, 0.025, 0.03.(b) shows plots of Im[F o m (ω)] using x = 0.0001, 0.001, 0.005, 0.01.The YSR bands overlap for x > 0.01 and Im[F o m (ω)] becomes the superposition of two semicircles.For this reason curves with higher values of x were omitted in (b).It was important to consider a small broadening η in order to cause the lines to overlap completely.

FIG. 4 .
FIG. 4. Im[F om (ω)] obtained using the same parameters of Fig.2and η = 10 −3 , but with different values of x, x = 0.01, where the two Shiba bands do not overlap, and x = 0.03, where the two bands do overlap, giving rise to a resonance around ω = 0. Solid blue lines were obtained using the DMFT equations (3), while dashed red lines correspond to Eq. (8).