Nematic Fluctuations in an Orbital Selective Superconductor Fe$_{1+y}$Te$_{1-x}$Se$_{x}$

We present a systematic study of the nematic fluctuations in the iron chalcogenide superconductor Fe$_{1+y}$Te$_{1-x}$Se$_{x}$ ($0 \leq x \leq 0.53$) using the elastoresistivity technique. Near $x = 0$, in proximity to the double-stripe magnetic order of Fe$_{1+y}$Te, a diverging $B_{1g}$ nematic susceptibility is observed. Upon increasing $x$, despite the absence of magnetic order, the $B_{2g}$ nematic susceptibility increases and becomes dominant, closely following the strength of the $(\pi, \pi)$ spin fluctuations. Over a wide range of compositions ($0.17 \leq x \leq 0.53$), while the $B_{2g}$ nematic susceptibility follows a Curie temperature dependence (with zero Weiss temperature) at low temperatures, it shows deviations from Curie-Weiss behavior for temperatures higher than $50K$. This is the opposite of what is observed in typical iron pnictides, where Curie-Weiss deviations are seen at low temperatures. We attribute this unusual temperature dependence to a loss of coherence of the $d_{xy}$ orbital, which is supported by our theoretical calculations. Our results highlight the importance of orbital differentiation on the nematic properties of iron-based materials.

The intricate interplay between magnetism and nematicity in different families of iron-based superconductors has attracted great interest in the past few years [1][2][3]. In iron pnictides, magnetism and nematicity are tightly coupled; the antiferromagnetic transition is always coincidental with, or closely preceded by, a tetragonal-to-orthorhombic structural transition. The proximity of the two transitions can be naturally explained within the spin-nematic scenario, where the structural transition is driven by a vestigial nematic order arising from fluctuations associated with the antiferromagnetic stripe transition (see Fig. 1(b)) [4][5][6]. In iron chalcogenides, the coupling between magnetism and nematicity is less obvious. FeSe undergoes a nematic phase transition without any long-range magnetic order [7,8], which has been interpreted as evidence that the nematic order in FeSe is of orbital origin [9]. Nevertheless, spin stripe fluctuations do develop below the nematic transition [10], and static stripe order can be induced by hydrostatic pressure [11,12].
While there are ongoing debates on the mechanism by which nematicity forms without static magnetism in FeSe [13][14][15][16][17], Fe 1+y Te 1−x Se x provides another platform to approach this problem. As selenium is replaced by tellurium (i.e. x is changed from 1 to 0), the nematic phase transition is suppressed [18], and inelastic neutron scattering experiments revealed a complex evolution of the spin correlations associated with different magnetic patterns [19][20][21][22]. In particular, close to optimal doping (x ∼ 0.5), the wave-vector of spin fluctuations at low temperatures is (π, π) [in the crystallographic Brillouin zone], identical to the antiferromagnetic order in the iron pnictides. As the tellurium concentration increases, both superconductivity and the (π, π) spin fluctuations disap- pear. The latter are replaced by short-range magnetic correlations near (0, π) that eventually condense into the static double-stripe phase in Fe 1+y Te [23] (Fig. 1(a)). Previous elastoresistivity measurements revealed a diverging B 2g nematic susceptibility in optimally doped Fe 1+y Te 1−x Se x , consistent with the existence of (π, π)  spin fluctuations [24]. This finding suggests that nematic and magnetic fluctuations remain strongly intertwined even in the absence of static nematic and magnetic orders. Nevertheless, in contrast to the magnetic sector, the behavior of nematic fluctuations for doping concentrations beyond optimal is still poorly characterized. The compositional dependence of the nematic susceptibility in Fe 1+y Te 1−x Se x would therefore constitute an important step in the effort to elucidate the relationship between nematicity and magnetism. Another motivation to study Fe 1+y Te 1−x Se x is to understand the influence of orbital selectivity on the nematic instability. Orbital selectivity (or orbital differentiation) refers to the fact that different orbitals are renormalized differently by electronic correlations, a characteristic property of Hund's metals that appears to be much more prominent in the iron chalcogenides in comparison with the pnictides [25][26][27]. Experimentally, recent scanning tunneling microscopy (STM) measurements revealed the impact of orbital differentiation on the superconducting state [28]. Theoretically, it has been suggested that orbitally selective spin fluctuations may be the origin of nematicity without magnetism in FeSe [29]. Nematic order was also proposed to enhance orbital selectivity by breaking the orbital degeneracy, leading to asymmetric effective masses in different dorbitals [30]. The effect of orbital differentiation becomes even more extreme as selenium is replaced by tellurium. In Fe 1+y Te 1−x Se x , angle resolved photoemission spectroscopy (ARPES) revealed a strong loss of spectral weight of the d xy orbital at high temperatures, which was interpreted in terms of proximity to an orbital-selective Mott transition [31]. Similar drastic changes were also observed as a function of doping [32], mimicking the evolution of spin fluctuations. Nevertheless, the impact of orbital incoherence on nematicity remains little explored [33].
In this report, we present systematic measurements of both the B 1g and B 2g nematic susceptibilities of Fe 1+y Te 1−x Se x (0 ≤ x ≤ 0.53) using the elastoresistivity technique. We demonstrate that the doping dependence of the two nematic susceptibilities closely track the evolution of the corresponding spin fluctuations. In particular, a diverging B 1g nematic susceptibility is observed in the parent compound Fe 1+y Te, suggesting that the spinnematic paradigm also applies to the double-stripe AFM order [34][35][36]. A diverging B 2g nematic susceptibility is observed over a wide range of doping (0.17 ≤ x ≤ 0.53), and its magnitude is strongly enhanced by both Se doping and annealing. In addition, the temperature dependence of the B 2g nematic susceptibility shows significant deviation from Curie-Weiss behavior above 50K. This is in sharp contrast to the iron pnictides, where the Curie-Weiss temperature dependence extends all the way to 200K. This unusual temperature dependence is captured by a theoretical calculation that includes the loss of spectral weight of the d xy orbital, revealing its importance for B 2g nematic instability.
Single crystals of Fe 1+y Te 1−x Se x were grown by the modified Bridgeman method. The electrical, magnetic and superconducting properties of Fe 1+y Te 1−x Se x are known to sensitively depend on y, the amount of excess iron. To study these effects, crystals were annealed in selenium vapor to reduce the amount of excess iron. By symmetry, the B 1g and B 2g nematic susceptibilities are proportional to the elastoresistivity coefficients m 11 −m 12 and 2m 66 , respectively. We performed the elastoresistivity measurements in the Montgomery geometry, which enables simultaneous determination of the full resistivity tensor, hence the precise decomposition into different symmetry channels, as illustrated in Fig. 1(c) and (d). Details of the Montgomery elastoresistivity measurements can be found elsewhere [24]. The crystal orientation is determined by polarization resolved Raman spectroscopy. Representative data of anisotropic resistivity as a function of anisotropic strain at 20K in B 1g and B 2g channels are shown in Fig. 1 66 shows a strong temperature dependence that grows continuously as temperature decreases. For x = 0.45, 2m 66 reaches a value of ∼ 100, comparable to optimally doped pnictides. While preserving a similar diverging temperature dependence, the maximum absolute value of 2m 66 decreases rapidly as selenium concentration decreases, from 100 for x = 0.45 to 8 for x = 0.17. On the other hand, m 11 − m 12 shows a diverging response when x is below 0.17, which is in the vicinity of the double-stripe AFM order. As selenium concentration increases, m 11 − m 12 evolves to a temperature independent response, with small kinks at low temperatures likely coming from contamination of 2m 66 due to misalignment. Overall, our observation of the doping dependence of 2m 66 and m 11 − m 12 is consistent with the evolution of low-temperature spin fluctuations from predominantly (π, 0) at small x to predominantly (π, π) at optimal doping x ∼ 0.5 [19][20][21][22].
To gain more insight, we fit the 2m 66 and m 11 − m 12 to a Curie-Weiss temperature dependence: For the parent compound Fe 1+y Te, m 11 −m 12 can be well fitted to a Curie-Weiss behavior in the temperature range just above the double stripe AFM ordering temperature T mag = 71.5K (Fig. 2(c)). The fitted Curie-Weiss temperature T * is slightly smaller than T mag . Despite the smaller absolute value (∼ 10 at maximum), the behavior of m 11 − m 12 is reminiscent of the 2m 66 in the parent phase of iron pnictides, suggesting that the spin-nematic mechanism is still at play here, in agreement with theoretical expectations [34][35][36]. temperatures, as can be seen in the linear temperature dependence of |2m 66 −2m 0 66 | −1 below 50K. It shows a significant deviation for temperature greater than 50K. The T * obtained from the low-temperature fitting is close to 0K. Intriguingly, the T * extracted from the Curie-Weiss fitting is approximately zero for all 0.17 ≤ x ≤ 0.53, while the Curie constant λ/a decreases with x (SOM). While the number of doping concentrations studied in the current work is insufficient to support a power law analysis, 2m 66 at constant T = 16K appears to be diverging as x increases from 0.17 to 0.45 ( Fig. 3(a)). Both the doping dependence and the near zero T * are consistent with the existence of a putative nematic quantum critical point at x ∼ 0.5. Interestingly, recent work doping FeSe with Te suggests that the 90K nematic transition of FeSe is continuously suppressed and extrapolates to 0 at x ∼ 0.5 [18].
This deviation from Curie-Weiss at high temperatures is very unusual. In the iron pnicitides, such a deviation was only observed at low temperatures in transitionmetal doped BaFe 2 As 2 ( Fig. 2(e)) and LaFeAsO. This unusual temperature dependence of 2m 66 appears to echo the coherent-incoherent crossover observed by ARPES [31], where the spectral weight of the d xy orbital is strongly suppressed as the temperature increases or as the selenium concentration decreases. To further confirm this correlation, we measured the Hall coefficient Inset of (a) shows the temperature dependence of the zero field cooling (ZFC) magnetic susceptibility measured at 100Oe (H ab). The superconducting volume fraction is significantly enhanced for the annealed sample.
R H , which has been demonstrated to be a good indicator of this incoherent-to-coherent crossover [32,37,38]. The recovery of the d xy spectral weight is generally correlated with a sign-change of R H [38] from positive to negative. Fig. 3(a) shows the low-temperature R H and 2m 66 as a function of doping, whereas Fig. 3(b-d) contain the full temperature and doping dependence of R H , m 11 − m 12 , and |2m 66 |, respectively. These plots reveal the strong correlation between a negative R H and an enhancement of 2m 66 . The properties of Fe 1+y Te 1−x Se x also depend on the amount of excess iron, which can only be removed by annealing [39]. Taking x = 0.45 as an example, the resistance of the annealed sample is metallic for temperatures below 150K (Fig. 4(a)). As Fig. 4(b) shows, at around 40K the Hall coefficient of the annealed sample turns from positive to negative, which is a signature of incoherent to coherent crossover. In contrast, the resistance of the as-grown sample shows a weakly insulating upturn at low temperatures ( Fig. 4(a) black dashed curve), and the Hall coefficient remains positive at all temperatures ( Fig. 4(b) black circles), indicating that the d xy orbital is still incoherent at low temperatures. Interestingly, at the same temperature where the resistance and the Hall coefficient of the as-grown and annealed samples depart from each other, the elastoresistivity coefficient 2m 66 shows a pronounced enhancement for the annealed sample ( Fig.  4(c)). Such an enhancement was observed in all annealed samples (SOM), providing further evidence of the correlation between the enhancement of the nematic susceptibility and the coherence of the d xy orbital. The doping and annealing dependences of 2m 66 presented above strongly suggest that the B 2g nematic susceptibility also have an orbitally-selective character. Indeed, previous theoretical works have highlighted the impact of orbital degrees of freedom on spin-driven nematicity [16,29,30,[40][41][42]. Using a slave-spin approach, Ref. [33] found a suppression of the orbitalnematic susceptibility due to orbital incoherence. To model our data, we employ the generalized random phase approximation (RPA) of Ref. [43] to compute the spin-driven nematic susceptibility for the five-orbital Hubbard-Kanamori model (details in the SM). For fully coherent orbitals, it was found that the largest contribution to the nematic susceptibility χ nem comes from the d xy orbital. Thus, one expects that χ nem would be suppressed if the d xy orbital were to become less coherent.
To verify this scenario, we calculated how χ nem changes upon suppressing the spectral weight Z xy of the d xy orbital. For our purposes, the reduction in Z xy acts phenomenologically as a proxy of the incoherence of this orbital, similarly to [28], but its microscopic origin is not important. Fig. 5(c)-(d) contrasts the nematic susceptibility for 0.7 ≤ Z xy ≤ 1. We note two main trends arising from the suppression of d xy spectral weight: first, as anticipated, the nematic susceptibility (and the underlying nematic transition temperature, which is non-zero in the model) is reduced (Fig. 5(c)). Second, its temperature dependence changes from a Curie-Weiss-like behavior over an extended temperature range to a behavior in which the inverse nematic susceptibility quickly saturates and strongly deviates from a linear-in-T dependence already quite close to the nematic transition ( Fig. 5(d)). These behaviors are remarkably similar to those displayed by the elastoresistance data shown in Fig.  5(a)-(b), with Z xy = 1 mimicking the behavior of optimally P-doped BaFe 2 As 2 and Z xy < 1, of optimally doped Fe 1+y Te 1−x Se x . Interestingly, the susceptibility associated with (π, π) fluctuations is also suppressed by the decrease in Z xy , in qualitative agreement with the neutron scattering experiments [44] (for a more detailed discussion, see SM). Of course, since Z xy in our model is an input, and not calculated microscopically, our model is useful to capture tendencies, but not to extract the experimental value of Z xy . Furthermore, note that in our calculation Z xy is temperature-independent, while in the experiment it changes with temperature.
In summary, our results reveal the close connection between nematic fluctuations and spin fluctuations in Fe 1+y Te 1−x Se x for both B 1g and B 2g channels. Additionally, the unusual temperature dependence of the B 2g nematic susceptibility can be attributed to the coherentto-incoherent crossover experienced by the d xy orbital, providing direct evidence for the orbital selectivity of the nematic instability. Our work presents Fe 1+y Te 1−x Se x as an ideal platform to study the physics of intertwined orders in a strongly correlated Hund's metal.
We thank Ming Yi for fruitful discussion.  warrants future investigation.

S3 Elastoresistivity measurement
The elastoresistivity measurement is conducted by gluing the square shape samples on a piezoelectric stack with the Montgomery contact geometry. Square edge of the samples are cut along the Fe-Fe (Fe-Ch) bonding direction for B 2g (B 1g ) elastoresistivity measurements, B 2g , respectively. One may also notice nonlinear responses in the elastoresistivity of annealed Fe 1+y Te 0.55 Se 0.45 are present in both the B 2g and A 1g channels. According to the symmetry constraints, the nonlinear term in the A 1g channel is most likely due to the second-order B 2g strain: While for the B 2g response (FIG. S3 (b)), the nonlinearity in the B 2g channel can be caused by either a mixed-in isotropic resistivity due to the sample's deviation from a perfect square along the Fe-Fe bonding direction, or a second-order strain B 2g A 1g .
Here, ( ∆ρ ρ 0 ) m B 2g is the measured value of resistivity anisotropy, which contains the mixed-in isotropic component. In this measurement, we mainly focus on the anisotropy term m B 2g B 2g , also known as 2m 66 in the Voigt notation, which is proportional to the nematic susceptibility in the B 2g symmetry channel of D 4h symmetry group. By finding the anisotropic strain neutral point and fitting the anisotropic elastoresistivity response quadratically, we are able to remove the mixed-in error and extract the first order elastoresistivity coefficient 2m 66 .
Here µ and ν are orbital indices, σ labels spin and we assume U = U − 2J and J = J. µν (k) denotes the dispersion, in this case obtained from a tight-binding fit to DFT. In the results shown in the main text, we used the band structure parameters presented in Ref. [2], which give three hole-like Fermi pockets and two electron-like Fermi pockets. While this tight-binding parametrization is not intended to model specifically optimally doped FeTe 1−x Se x , it offers a solid framework to elucidate, on general grounds, the tendencies of how the nematic susceptibility is affected by the suppression of d xy orbital spectral weight.
The interactions can be conveniently expressed as elements of the same tensor: where repeated indices are not summed.
To elucidate how orbital differentiation affects the nematic susceptibility within our framework, we follow the approach outlined in Ref. [6] and modify the Green's function by: where we stress that, once again, the Einstein convention is not assumed. Since both g µνρλ and χ µν (q) depend on the Green's functions, a reduced Z factor will have an intricate impact on the nematic susceptibility. As discussed in the main text, this is a phenomenological way to mimic the complicated effect of incoherence and loss of spectral weight on the response function, whose validity requires that the system remains in a metallic state.
The nematic susceptibility plotted in the main text for different values of Z xy corresponds to the largest eigenvalue λ nem of the equation: In Fig. S5, we show how the corresponding eigenvector Φ (n) µν changes for decreasing Z xy . Here we used U = 1.2 eV and J = U/6. To highlight the impact of Z xy we fixed the filling to 5.9 electrons per site. As expected, the main contribution to the nematic susceptibility, signaled by the brightest squares in the figure, shifts from intra xy-orbital processes to intra xz/yz-orbital processes.
A direct consequence of the reduction of the nematic susceptibility is a suppression of the nematic transition temperature, which in our model is manifested as a divergence of the largest eigenvalue λ nem . In Fig. S6, we plot both the nematic transition temperature T nem and the bare (i.e. non-renormalized by nematic order) magnetic transition temperature T mag as a function of Z xy . Note that not only are both transition temperatures strongly suppressed, but their separation also decreases significantly for decreasing Z xy . As noted in Ref. [3], these RPA transition temperatures are, not surprisingly, overestimated with respect to the actual transition temperatures. For this reason, and to be able to compare the temperature dependencies of the nematic susceptibilities of systems with very different values of T nem , in the main text we plot λ nem as a function of T − T nem .
The suppression of T mag for decreasing Z xy is also manifested in the suppression of the overall magnitude of the bare magnetic susceptibility calculated at the spin-stripe wavevector, χ mag , as shown in Fig. S7. Here, χ mag was calculated in a manner similar to the µν for Z xy = 1 and Z xy = 0.7 at a temperature immediately prior to the nematic instability of each case. For Z xy = 1, the dominant contribution to the nematic susceptibility arises from the xy orbital. As Z xy is reduced, the xz and yz orbitals become the dominant ones. Nematic (blue) and bare magnetic (red) transition temperatures as a function of 1 − Z xy . Both are reduced by reducing Z xy , along with their relative separation. At Z xy = 0.7 the separation between the two vanishes within our temperature resolution (δT < 0.2 meV).
nematic susceptibility: where, in Fig. S7, q = Q 1 or q = Q 2 and we show the leading eigenvalue. Note that, in contrast to the nematic susceptibility, the temperature dependence of χ mag is not strongly affected by the reduction of d xy spectral weight, as seen clearly from Fig. S7. However, for a fixed temperature, there is a strong suppression of χ mag with decreasing Z xy . This last feature is in qualitative agreement with neutron scattering experiments, which showed that, in optimally-doped FeTe 1−x Se x , the low-energy stripe-type magnetic fluctuations are suppressed with increasing temperature [7]. According to our analysis in the main text, increasing the temperature promotes a less coherent d xy orbital. Note, however, that Ref.
[7] also reported that, as temperature was increased, besides a suppression of intensity at the stripe wave-vector, an incommensurate peak appeared in the momentum-resolved magnetic susceptibility at low energies. We did not observe such incommensurate peaks in our energy integrated magnetic susceptibility, i.e. the peak remains at q = Q 1 (and at q = Q 2 ), suggesting that this effect cannot be captured phenomenologically by a constant Z xy factor.