Observation of orbital ordering and origin of the nematic order in FeSe

To elucidate the origin of nematic order in FeSe, we performed field-dependent 77Se-NMR measurements on single crystals of FeSe. We observed orbital ordering from the splitting of the NMR spectra and Knight shift and a suppression of it with magnetic field B0 up to 16 T applied parallel to the Fe-planes. There is a significant change in the distribution and magnitude of the internal magnetic field across the orbital ordering temperature Torb while stripe-type antiferromagnetism is absent. Giant antiferromagnetic (AFM) spin fluctuations measured by the NMR spin-lattice relaxation are gradually developed starting at ~ 40 K, which is far below the nematic ordering temperature Tnem. These results demonstrate that orbital ordering is the origin of the nematic order, and the AFM spin fluctuation is the driving mechanism of superconductivity in FeSe under the presence of the nematic order.


Introduction
The interactions between structure, magnetism and superconductivity in Fe-based superconductors have been of wide interests [1][2][3]. The experimental determination of these interactions is challenging due to the occurrence of nematic order often at or near the temperature of a stripe-type long-range antiferromagnetic (AFM) order [4][5][6][7][8]. Similar to the stripe-type AFM order, the nematic order also breaks the lattice four-fold (C 4 ) rotational symmetry of a high-temperature phase, as evidenced by a tetragonal-to-orthorhombic structural phase transition at T s [7][8][9][10][11]. On the other hand, the nematic order is directly linked to the superconducting state because nematic instability is a characteristic feature of the normal state upon which at lower temperatures the superconductivity emerges [1,8,12,13], and thus nematicity is deemed a precursor of superconductivity in unconventional superconductors including the cuprates. It is generally believed that the nematic order is electronic and the structural phase transition is the consequence of the nematic order, since the lattice distortion is much smaller than the observed anisotropy of the in-plane resistivity in the nematic phase [9,14]. However, it remains highly controversial regarding the origin of the nematic order whether it is driven by spin order [15,16], AFM spin fluctuations [15][16][17], and/or orbital order [18][19][20][21][22][23][24].
However, recent findings show that stripe-type AFM order emerges under high pressure, and the AFM ordering temperature increases with high pressure [17,33,[35][36][37]. These findings make the origin of the electronic nematic order more elusive. Even though various experimental techniques have been used for the study, most research work reported was focused on the doping and high pressure effects on the properties of FeSe. A systematical investigation of the effect of applied magnetic field on the properties of FeSe is still lacking.
Here we present for the first time field-dependent 77 Se-NMR measurements systematically on high-quality single crystals of FeSe with applied magnetic field B 0 up to 16 T and temperature down to 1.5 K. We observed orbital ordering which is demonstrated by the splitting of the NMR spectrum and Knight shift. As measured by the 77 Se-NMR spin-lattice relaxation, giant AFM spin fluctuations gradually develop starting at ∼ 40 K, which is far below the nematic order temperature T nem = T s = 89 K. These discoveries provide direct evidence of orbitaldriven nematic order in FeSe. They also shed light on the important role of the nematic order on the superconductivity of Fe-based superconductors.

Experimental Section
Single crystals of FeSe were grown in evacuated quartz ampoules using the AlCl 3 /KCl flux technique with a temperature gradient of 5 o C/cm along the ampoule length. The temperatures of the hot and cold ends used for the growth were 427 o C and 380 o C, respectively. X-ray diffraction verified that the crystals have a high-purity single phase with a tetragonal crystal structure at room temperature, where the lattice c-axis is perpendicular to the Fe-planes (ab-plane). A SQUID magnetometer was used to measure the DC magnetic susceptibility χ(T ). The samples used for our NMR measurements have a typical size of 3.3 × 2.7 × 0.1 mm 3 .
The NMR coil used for the measurements was made from 50 µm diameter silver wire wound with ∼ 18 turns and attached to a goniometer on an NMR probe by epoxy. A single-crystal FeSe sample was put inside the coil so that the sample rotation axis is in the lattice abplane and perpendicular to the applied field B 0 . A commercial network analyzer was used for the observation of the tuning and matching of the resonant circuit located at the bottom of the NMR probe. The NMR spectra and spin-lattice relaxation time T 1 were measured with an inversion-recovery method, where a π pulse is first applied to invert the nuclear magnetization M 0 to the -z axis, and then after a delay time t, a π/2 pulse is applied to measure the recovering magnetization M (t) component along the z axis, which gives T 1 as a function of time t as M (t)

Results and Discussion
Figure 1(a) shows the typical 77 Se-NMR spectra at B 0 = 12 T and temperature T = 40 K (below T nem ), by the variation of the angle θ between B 0 and the lattice c-axis of FeSe. The spectrum splits into two peaks (P 1 and P 2 ) which is observed when B 0 is applied ∼ 25 o from the c-axis. The splits reach the maximum when B 0 a&b, indicating the largest anisotropy of the internal field at the Se-sites in the a&b-plane.  Noticeably, the 77 Se-NMR spectra are fully magnetic with no electron charge or quadrupolar contributions because 77 Se is a spin I = 1/2 nucleus (which has no quadrupole moment). The spectral splitting was only observed when temperature T is less than T nem = T s = 89 K, but not for B 0 c (the reason to be revealed later).
Thus, undoubtedly the spectrum split is the result of a structure symmetry break in the ab-plane due to the tetragonal-to-orthorhombic structure phase transition (structurally a and b are not equal any more), which is known as the consequence of the electronic nematic order in the Fe-planes [23,24].
Moreover, with the nematic order the spectrum splits also reflect a significant change in spacial field distribution (∆B FWHM ) and also a change in the value of the internal field (B ′ ) at the Se-sites. Here ∆B FWHM = ∆f / 77 γ I , where 77 γ I = 8.131 MHz/T is the gyromagnetic ratio of the 77 Se nucleus, and B ′ = (ν − ν L )/ 77 γ I , where ν is the NMR frequency and ν L is the Larmor frequency (ν L = 77 γ I B 0 ). the Se-sites has a change ∆B ′ = ± 12.0 G (a Knight shift change of ± 0.010%) from the average value (B ′ ) of the internal field B ′ = 160 G (or an average Knight shift 0.133%) (not shown here), i.e., the change of the value of internal field ∆B ′ reaches ± 7.5% from the average value of the internal field B ′ in the Fe-planes.
The Knight shift K is defined by K = (ν − ν L )/ν L as a tradition, and it should be field independent. That the values of K(T ) at B 0 a&b are apparently larger than those at B 0 c at T < T nem indicates an anisotropic hyperfine coupling.
In general, the Knight shift K is given by [39,40]: Here χ spin and χ orb are the electron spin and orbital susceptibility, respectively. A spin and A orb are the hyperfine coupling constants between the studied nucleus and the electron spins and the electron orbitals, respectively. N A is the Avogadro's number and µ B is the Bohr magneton. Likewise, the magnetic susceptibility χ is the sum of the contributions from core diamagnetic susceptibility (χ dia ), orbital (van Vleck) paramagnetic susceptibility (χ orb ) and Pauli spin paramagnetic susceptibility (χ spin ) [41,42], i.e., χ = χ dia + χ orb + χ spin . Here, for FeSe, χ dia = − 6.1 × 10 −5 cm 3 /mol from the diamagnetism of the atomic ions, and χ orb is T -independent unless there is an orbital change such as an orbital ordering. Interestingly, at T < T nem , K(T ) versus χ(T ) gradually deviates from the high temperature linear relation, as seen in figure 2(a) for both B 0 a&b and B 0 c. Because K(T ) and χ(T ) are fully magnetic in nature, as described by equation (1), this deviation can only be explained by a change in the electron spin susceptibility χ spin (T ) such as that as a result of an AFM order of the electron spins or AFM spin fluctuations, and/or by a change in the electron orbital susceptibility χ orb such as that as a result of an ordering of the electron orbitals, as well as associated changes in the hyperfine couplings to the electron spins (A spin ) and/or to the electron orbitals (A orb ), any of which could lead to a change in K(T ) simultaneously. This is also seen by the expression [39,40,42] where only K(T ) and χ(T ) are temperature dependent.
However, surprisingly, upon further cooling in temperature, the K(T ) -χ(T ) plot exhibited in figure 2(a) inset shows that the slope of K(T ) versus χ(T ) is ∼ 0, both at B 0 a&b and B 0 c at B 0 = 12 T in the temperature range ∼ 60 -10 K, which is a wide range of temperature below T nem and above the critical temperature T c of superconductivity, i.e., K spin, a ≈ 0, K spin, b ≈ 0, and K spin, c ≈ 0. This is also true for all other fields we applied. That is to say that below T nem the spin Knight shift K spin (T ) becomes negligible at all directions, i.e., K ≈ K orb .
In other words, the Knight shift K(T ) at low temperatures predominantly comes from the contribution of the orbital Knight shift K orb [ Figure 2 The reason that K spin (T ) ≈ 0 in all directions can be understood by enormous AFM spin fluctuations developed in the same temperature regime, whereas there is no existence of electron spin order, as directly evidenced by our 77 Se-NMR spin-lattice relaxation data (see next), with the consideration of a more general expression of the spin Knight shift as [39,42]: . This is the summation of the terms of hyperfine coupling interaction to the individual electron spins (the degree of electron spin polarization is ∝ χ i spin ), where each term could be very different from each other due to the AFM spin fluctuations, resulting in a cancelation of them.
On the other hand, the dramatic increase of the orbital Knight shift K orb [ Figure 2(b)] below T S (T nem ) must be the result of an orbital ordering. To confirm this, we studied the internal field difference (∆B ′ a,b ) in the ab-plane by the measurement of the frequency difference (∆ν a,b ) of the NMR spectrum peaks (P 1 and P 2 ), as shown in figure 2(c). ∆ν a,b reaches ∼ 12.5 kHz and 25.0 kHz, or a value of internal field difference ∆B ′ a,b ≈ 15.6 G and 31.2 G at B 0 = 8 T and 16 T, respectively, at low temperatures. These values are scalable with B 0 , which is understandable as they are magnetic in nature. Here we have ∆B ′ a,b = ∆ν a,b / 77 γ Se . Since there is no appearance of AFM spin order and the in-plane anisotropy of the paramagnetic spin Knight shift is expected to be negligible (i.e., K spin, a ≈ K spin, b ), from the Knight shift we have (3) Therefore, this indicates that all the data values of ∆ν a,b shown in figure 2(c), are essentially completely from the orbital contributions (for convenience, we say all orbital), i.e., the internal field difference in the abplane is fully determined by the difference of the hyperfine coupling to the Fe-electron orbitals among the a− and b−axes. In other words, these data verify that there is an electron orbital ordering occurring at T ≤ T nem . Now we can define a characteristic orbital ordering temperature T + by the intersection of two lines that fit to the data in the transition area as shown in figure 2(c), and we find that T + is linear to B 0 as: Figure 4]. We noted that T + indicates the temperature where orbital ordering is fully developed while the orbital ordering starts at T nem upon cooling. Here ∆ν a,b or ∆B ′ a,b can be treated as the orbital ordering parameter [23], ∆ν a,b ∝ √ T nem − T near T = T nem , and as B 0 → 0, T + = T nem .
Furthermore, in order to investigate the electron spin dynamics and to support the observations in the NMR spectrum and Knight shift, we performed the 77 Se-NMR spin-lattice relaxation measurements as a function of temperature and applied field, as exhibited in figure 3.
Generally, 1/T 1 T probes the imaginary part of the lowfrequency (ω → 0) dynamical susceptibility [χ(q, ω)] averaged over the momentum (q) space as [39,43] , which gives the Korringa law [39,43]: where γ I(e) is the gyromagnetic ratio of nucleus (electron), N (E F ) is the density of states of electrons at the Fermi energy E F , and f (E) is the energy distribution function. For AFM correlated electrons, χ(q) can have a peak at the AFM wave factor Q = (π, π), and then 1/T 1 T ∝ χ(Q) with a Curie-Weiss type relation as: 1/T 1 T = C ′ /(T − θ), as often seen in cuprate and other Fe-based superconductors [6,[44][45][46]. For AFM fluctuations, the fit parameter θ < 0, and for large spin fluctuations C ′ is large. Thus, important information can be obtained from the NMR spin-lattice relaxation. First, figures 3(a) and 3(b) show the nematic order/structure phase transition at T nem = T s , which is independent of B 0 . Second, enormous AFM spin fluctuations are developed [a significant deviation from the Korringa law] as seen in the plot of 1/(T 1 T ) versus T but they start at ∼ 40 K and below only [inset of figures 3(a) and 3(b)], which is far below T nem = 89 K. With the fit to the Curie-Weiss relation for 10 K < T < 40 K, we have the values of θ = −4.6 (−21.5) K, and C ′ = 10.0 (7.2) s −1 for B 0 a&b (B 0 c). Here θ is comparable while C ′ is much larger than those of other Fe-based superconductors [47][48][49][50]. Third, the AFM spin fluctuations drop significantly at T < T c due to diamagnetism associated with the pairing symmetry of the electron spins, and 1/T 1 ∝ T α , where α ≈ 3 in low fields, consistent with a line-node gap behavior of a d-wave superconductor, agreeing with reports on various Fe-based superconductors [35,40,47,51]. Moreover, below T nem in the range 40 K<T ≤T nem , 1/T 1 T continues to show a free-electron behavior (Korringa law), where 1/T 1 T is a constant as exhibited in figures 3(a)-3(b) insets, i.e., there is essentially no AFM spin fluctuations at the T nem regime over a wide range of temperature. Therefore, AFM spin fluctuations can be excluded from the driving mechanism of the nematic order. Figure 4 shows temperature -field phase diagram we obtained. First, as described in detail above, the applied magnetic field decreases the characteristic orbital ordering temperature T + rather sensitively. Second, the structural phase transition temperature T s and the nematic ordering temperature T nem are not affected by the applied field (T nem = T s = 89 K). Third, stripe-type AFM order is absent at all the applied fields.
Here the values of T c were determined by our resonance frequency method [54]. We would like to point out that, among the two groups of the 3d Fe t 2g (d xy , d yz , d xz ) and e g (d x 2 −y 2 , d 3z 2 −r 2 ) orbitals (totally five orbitals), d xy , d x 2 −y 2 , and d 3z 2 −r 2 are rotationally symmetric in the xy-plane. Thus, the only two candidates related to the tetragonalorthorhombic degeneracy (lattice symmetry breaking) are the d yz and d xz orbitals. And, what an NMR spectrum directly measures is the local field distribution and local field values parallel to the externally applied magnetic field at the nucleus at the atomic scale. Therefore, because of the unique symmetry of the the d yz and d xz orbitals in the lattice axis c-direction the ordering of these d yz and d xz orbitals is not able to cause any splitting of the 77 Se-NMR spectra at B 0 c in FeSe.
Finally, we discuss the field effect on the characteristic temperatures. That the values of T s (T nem ) are not affected by the directions or magnitude of the applied field could be explained by the weak anisotropy character of the paramagnetic Fe-spins in the high symmetry tetragonal lattice. That T + is linearly proportional to B 0 could be understood due to its full magnetic character that involves electron orbital moments, while the reason for the decrease of the value of T + with B 0 is not clear. Whereas there is no appearance of long-range AFM order of the electron spins, we note that there were reports about short-range stripe magnetic order and possible spin-orbital coupling [55][56][57] developing in the nematic phase. From our data of the temperature-and fielddependence of 77 Se-NMR spin-lattice relaxation, 1/T 1 T versus T , we clearly see that giant AFM fluctuations gradually develop starting at ∼ 40 K and below, which is far below the nematic order temperature T nem . And above T nem the data of 1/T 1 T well follows the Korringa law, which is a direct evidence of trivial electron correlations. These results leave the orbital ordering unequivocally as the dominant driving force of the nematic order.

Conclusions
In summary, we report direct observation of orbital ordering which is demonstrated by the splitting of the 77 Se-NMR spectrum and Knight shift in single crystals of iron-based superconductor FeSe. As illustrated by the field-dependence of the 77 Se-NMR spin-lattice relaxation, stripe-type AFM order is absent, whereas giant AFM spin fluctuations developed starting at temperatures far below the electronic nematic order temperature T nem , thus both of which can be unambiguously excluded from the origin of the nematic order. These discoveries provide direct evidence of orbital-driven nematic order in FeSe. Our results also help to the understanding of the strong interplay between structure, magnetism and superconductivity in Fe-based superconductors as well as other unconventional superconductors.