Multiple fermion scattering in the weakly coupled spin-chain compound YbAlO3

The Heisenberg antiferromagnetic spin-1/2 chain, originally introduced almost a century ago, is one of the best studied models in quantum mechanics due to its exact solution, but nevertheless it continues to present new discoveries. Its low-energy physics is described by the Tomonaga-Luttinger liquid of spinless fermions, similar to the conduction electrons in one-dimensional metals. In this work we investigate the Heisenberg spin-chain compound YbAlO3 and show that the weak interchain coupling causes Umklapp scattering between the left- and right-moving fermions and stabilizes an incommensurate spin-density wave order at q = 2kF under finite magnetic fields. These Umklapp processes open a route to multiple coherent scattering of fermions, which results in the formation of satellites at integer multiples of the incommensurate fundamental wavevector Q = nq. Our work provides surprising and profound insight into bandstructure control for emergent fermions in quantum materials, and shows how neutron diffraction can be applied to investigate the phenomenon of coherent multiple scattering in metals through the proxy of quantum magnetic systems.

T he one-dimensional (1D) antiferromagnetic (AFM) XXZ spin S = 1/2 chain in a magnetic field (Eq. (1)) is among the most attractive, simple and, therefore, well-studied models in solid-state physics, from both experimental [1][2][3][4] and theoretical [5][6][7] points of view. Its Hamiltonian reads where Δ parameterizes the exchange anisotropy and H z refers to the applied longitudinal magnetic field. Depending on H z , the XXZ model exhibits three types of low-lying excitations: (i) fractionalized fermionic spinons at zero field; (ii) bosonic magnons in the field-polarized regime; (iii) complex quasiparticles known as psinons and psinon-antipsinons at intermediate field, which adiabatically connect spinons and magnons, respectively 8 . Note that though the psinons are not canonical fermions, their low-energy dispersion can be linearized, thus the low-energy excitations can be considered as right and left-moving fermions. These fermions are fractionalized quasiparticles, which means that any local operator (such as a spin operator, responsible for neutron-sample interaction in the inelastic neutron scattering (INS) cross-section) creates more than one excitation, and therefore their single-particle bandstructure is not directly accessible by spectroscopic techniques. However, the multiparticle excitation spectrum manifests itself in the dynamical spin susceptibility, χ″(q, ω), which can be probed by INS.
Having a system with a gapless spectrum of fermionic quasiparticles, one can predict its physical properties using approaches developed for conventional metals. In particular, the realization of the Fermi surface (FS)-the most important concept in the physics of metals-was discussed widely in the context of quantum spin liquids [9][10][11][12] and Heisenberg or XY spin chains, and in the latter case, the FS is reduced to two points at ±k F , as represented in Fig. 1a-c. Compared with electrons in a metal, the fermions in a spin chain have a significant advantage for experimentalists-by applying a magnetic field one can control the chemical potential directly without the doping-induced disorder, the need to perform complicated gating, and other problems intrinsic to electronic systems. The magnetic field effectively shifts the chemical potential of the fermion band and the nesting wavevector connecting the two branches follows directly the field-induced magnetization, 2k F = π(1 ± 2m) (where m is the magnetization normalized by its saturation value, M/M S ) 13 . In the Ising regime (Δ > 1), the 2k F scattering gives rise to a gapless incommensurate longitudinal mode at finite fields 14 . This mechanism is similar to FS nesting and stabilizes an incommensurate spin-density wave (SDW) order at q = 2k F when the interchain coupling is present. Even though theoretical predictions for this SDW ordering have existed for a long time, its experimental realization is limited to very few materials including BaCo 2 V 2 O 8 and SrCo 2 V 2 O 8 15,16 , both having a moderate Ising anisotropy, Δ ≃ 2.
Recently, field-induced SDW order was observed in the Heisenberg spin-chain material YbAlO 3 17 . In this compound, the combination of strong spin-orbit coupling (SOC) and crystalline electric field lifts the degeneracy of the ground-state multiplet and causes an effective S = 1/2 at each Yb site at low temperatures. Despite the strong SOC of the Yb ion, the dominating intrachain exchange interaction J = 0.21 meV is almost isotropic, with Δ ≃ 1 17,18 and the excitation spectrum above T N exhibits spinon continuum as shown Wavevector along the chain direction Wavevector along (0,0,L) q = π-δ 2 q = π+δ 1 Orthogonal wavevector (0,K,0), (r.l.u.)  in Fig. 1d. The finite dipolar interchain coupling with an Ising-like anisotropy 19 causes a commensurate AFM order below T N = 0.88 K.
The field-temperature phase diagram is summarized in Fig. 2a.
Application of a longitudinal magnetic field (B∥a) suppresses the commensurate ordering at low temperatures in favor of an incommensurate (IC) phase [ Fig. 1e, f] at B c1 = 0.32 T, and drives a quantum phase transition to the field-polarized phase at B c2 = 1.15 T. Moreover, the magnetization curve shows a weak plateau close to m = 1/3, as shown in Fig. 2b. We note that the observed incommensurate SDW order in YbAlO 3 is surprising, because the lowenergy physics of a Heisenberg chain is described by a Tomonaga-Luttinger liquid (TLL) in which the 2k F scattering is irrelevant at finite fields. Although it has been recognized that the Ising anisotropy of the interchain coupling has a crucial role 20,21 , this fieldinduced SDW order in YbAlO 3 remains poorly studied.
In this work, we have performed elastic and INS measurements to reveal the details of magnetic ordering within the lowtemperature SDW phase of YbAlO 3 . We found that the combination of isotropic intra-and Ising-like interchain interactions gives rise to the consecutive formation of the SDW and transverse antiferromagnetic (TAF) phases as the longitudinal magnetic field is increased. The interchain coupling within the SDW produces an effective incommensurate potential, which opens a way to multiple scattering of the fermions (the analogy of Umklapp process in metals and graphene [22][23][24]. This results in the formation of satellites at multiples of the incommensurate fundamental wavevector, nq, which we have observed in the neutron diffraction experiments. To verify the physical picture of multiple scattering we have performed detailed numerical calculations by density matrix renormalization group (DMRG) on quasi-1D systems and by QMC on the 3D coupled-chain model.

Results and discussion
We begin the presentation of our results with the INS spectra measured in the IC phase using Cold Neutron Chopper Spectrometer (CNCS).   Red arrows indicate extrinsic Bragg reflections due to the presence of the twin in the sample, which is discussed in detail in Sec. (S1) of the SI. Subtraction of the 2 T data set was applied to all spectra in c1-f2.
with a magnetic field. Previously, it was shown that the essential features of the spectra obtained at B ≥ 0.6 T can be reproduced by a S = 1/2 Heisenberg spin-chain model in a magnetic field without introducing any kind of interchain coupling, J ab 17 . This is because the interchain coupling has a dipolar origin, and is therefore proportional to the ordered moment, which is suppressed with increasing field. Only the spectrum collected at B = 0.4 T shows qualitatively new features at high energies and the introduction of interchain coupling is essential to obtain its correct description, but J ab affects the low-energy part of the spectrum, which was used to extract 2k F only very weakly. Thus, by examining the spectra at ℏω ≈ 0 + ε (ε is a small number, on the order of the spectrometer resolution) we can gain information about the quasi-zero-energy response of 1D subsystem as shown in Fig. 1e, f. However, even though J ab is weak, it is an essential ingredient of the system, because it induces the 3D long-range order, which manifests itself by the appearance of sharp Bragg peaks in the elastic channel. Let us turn to the description of the field-induced magnetic order. As discussed in the introduction, the magnetic field effectively shifts the chemical potential of the fermion band and induces a nesting of the FS at q L = 2k F = π(1 ± 2m). In the spectrum, this instability corresponds to the zero-energy position of the longitudinal mode, S zz (2k F , ω = 0 + ε), which we call the soft mode throughout the text. Within the SDW phase at B = 0.4 and 0.6 T, the soft-mode positions, 2k F = 1 ± 0.13 and 1 ± 0.28 (r.l.u.), respectively, can be extracted from the inelastic spectrum (see details of the soft-mode fitting in Sec. (S3) of the Supplementary Information (SI)) and the elastic Bragg peak is located exactly at the same position, q L1 = 1 ± δ 1 . When the field is increased to 0.8 T, the zero-energy mode position, 2k F = 0.59 (r.l.u.), shifts continuously away from (001), but the position of the elastic satellite, q L1 = 0.66 (r.l.u.), deviates from 2k F . At B = 1 T, the first set of satellites disappears, whereas the zero-energy mode position shifts gradually towards the ferromagnetic (FM) (002) and (000) positions. We summarize the evolution of the soft mode as a function of the magnetic field in Fig. 3a. The backgroundsubtracted INS cuts along the (00L) direction, integrated just above the elastic line, are plotted at different magnetic fields so that one can see the almost linear evolution of the soft-mode position, 2k F .
Note that some part of the intensity remains almost unchanged at the (001) position up to 0.75 T and then exhibits a second splitting, δ 2 . The field dependence of δ 2 is similar to δ 1 . However, the second splitting, as well as the overall behavior of the (001) peak, is not intrinsic to the spin-chain physics of YbAlO 3 , but rather is caused by the small degree of twinning in the sample, with the second twin rotated by 90 ∘ with respect to the main crystalline axes and oriented with [010]∥B (see Sec. (S1) of the SI for a detailed analysis of the second splitting). Henceforth, we will   c Representative ENS scans along the (00L) direction measured on FLEXX at T = 0.05 K at different magnetic fields. Red arrows mark extrinsic reflections due to the twin. d S zz (q, ω = 0) within the spin-density wave (SDW) phase calculated using quantum Monte-Carlo (QMC). e S zz (q, ω = 0) calculated by density-matrix renormalization group (DMRG) method for different external magnetic fields with (left) and without (right) the incommensurate field, h ic . One can see that h ic stabilizes the ordering at the fundamental wavevector 2k F and produces additional harmonics at 2nk F . f1-f2 Distorted SDW with square-like modulation calculated by DMRG (f1) compared with a harmonic SDW (f2). The distorted SDW in f1 was calculated using Eq. To further reveal the details of the magnetic structure of the IC phase, we performed an elastic neutron scattering (ENS) experiment using the triple-axis spectrometer (TAS) FLEXX 25 and the primary results of our measurements are summarized in Fig. 3a, c. Figure 3c shows the raw ENS scans measured along the (00L) direction as a function of magnetic field at T = 0.05 K. One can see that below the first critical field, B c1 , the diffraction patterns contain commensurate AFM and structural reflections at (001) and (002), respectively. Above B c1 , the AFM peak splits into q L1 = 1 ± δ 1 , and δ 1 follows the position of the soft mode of the spectra and magnetization curve up to the 1/3 plateau. Surprisingly, we found that, in addition to the primary IC satellites at q L1 , the diffraction pattern exhibits two additional weak peaks, whose positions evolve with magnetic field as q L2 = 2 ± 2δ 1 (B) and q L3 = 1 ± 3δ 1 (B), as shown in Fig. 3c. The intensities of those peaks are between 100 and 1000 times weaker than the structural and primary IC peak. For this reason we were able to resolve only the strongest harmonics of q L2 at B = 0.4 T in our CNCS measurements, visible as the weak peak at (0 0 0.27) r.l.u. in Fig. 2c2.
At B = 0.7 T, the q L1 and q L2 peaks merge at ð0 0 4 3 Þ and the ENS signal can be deconvolved using two peak functions: (i) a sharp peak, whose position is locked at q L = 1.33 (r.l.u.) and whose intensity is suppressed gradually with field, until it finally vanishes above 0.85 T; (ii) a broad peak, which corresponds to the soft-mode contribution and shifts continuously with magnetic field. We summarize the field dependence of the IC peaks along with the soft-mode position in Fig. 3a. The locking of the Bragg peak position takes a place around the same field range where the M S /3 plateau was observed. In spite of the presence of the soft mode at finite wavevector, q = 2k F , in the field range between 0.9 T and B c2 , the sharp incommensurate peak is absent, indicating that the SDW phase is suppressed, although the sample remains in an ordered phase, as is evident from the sharp features in thermodynamic measurements 17 . This suggests that the magnetic ordering changes to the transverse channel, which supports a gapless Goldstone mode as a consequence of the broken in-plane U(1) symmetry. The soft mode in the longitudinal channel then becomes overdamped and fades out with increasing field.
Indeed, a TAF phase with nonzero transverse component of the spin structure factor, S xx (q = π), for H z ≳ 0.7 T is confirmed in both (QMC) and DMRG calculations 20,21 . The TAF phase is continuously vanishing at B c2 = 1.15 T. Unfortunately, owing to the presence of the small twin in our sample, we could not compare the experimental field-dependence of the peak at q = π quantitatively with the theoretical predictions, but qualitatively both follow the same trend and gradually vanish close to B c2 (see Sec. (S1) of the SI for an analysis of the (0 0 1) peak intensity).
The origin of the satellites can be understood as a competition of intrachain Heisenberg and interchain Ising-like interactions. The latter enhance longitudinal spin correlations and therefore support a collinear ordering with maximal 〈S z 〉. On the other hand, at the finite field, the Heisenberg term promotes ordering in the transverse components. The compromise between these energy scales can be achieved by the formation of an incommensurate SDW with a weak modulation at integer multiples of the fundamental propagation vector nq.
To better understand this phenomenon, we first map the Heisenberg spin chain in a magnetic field (Eq. (1) with Δ = 1) to a system of interacting fermions (full details can be found in Sec. (S4) of the SI). The low-energy physics is dominated by the leftand right-moving fermions, respectively, Ψ L at −k F and Ψ R at k F . Then, following the standard bosonization procedure (see Sec. (S4) of the SI and ref. 17 ), the physics of these fermions can be understood by the following Luttinger model where u and K are the spin-wave velocity and the Luttinger parameter, respectively. This model describes a TLL with dominant transverse correlations. In terms of fermions, it shows that in addition to the band renormalization effect, the intrachain Ising interaction is irrelevant at low energies. However, the Ising term of the interchain coupling, which is introduced as an effective contribution to the field term, produces scattering at q = 2k F between the left-and right-moving fermions [ Fig. 3b]. This Umklapp process adds a term J ab cosð2ϕÞ to Eq. (2), meaning that it is a relevant perturbation of the Luttinger model, and its consequence is the band-folding of the left-and right-moving fermions with the incommensurate wavevector 2k F shown in Fig. 3b. In 3D coupled chains, this instability of FS gives rise to an SDW order at the incommensurate wavevector q = 2k F and opens a small gap at the Fermi level. Moreover, the nth order perturbation in J ab create processes involving multiple coherent (Umklapp) scattering of fermions with incommensurate wavevectors nq = 2nk F , which accounts for the observed satellite peaks in the ENS measurements.
To further demonstrate that the bosonization analysis describes the physics for the satellites, we perform numerical simulations on microscopic spin Hamiltonians. We first use DMRG to study an effective single-chain spin model, In this equation, the second term shows the effect of the external magnetic field while h ic denotes the interchain incommensurate "molecular" field produced by the neighboring chains; 2k F is determined for each H z as 2k F = π(1 − 2m). S zz (q, ω = 0) calculated for different H z in the presence and absence of modulated fields is shown in Fig. 3e. It is clear that, in absence of the incommensurate modulation, S zz shows only a broad feature at the fundamental wavevector, q = 2k F , whereas h ic = 0.4 J produces a series of additional satellites at nq. The resulting distortion of the SDW is represented in Fig. 3f1 has a pronounced "square"-like shape when compared with the harmonic SDW [ Fig. 3f2]. We note that the strength of internal field h ic used in DMRG calculations was overemphasized, compared with its quantitative value in YbAlO 3 , in order to present the type of distortion more clearly. Finally, to show that the molecular-field approximation used in our DMRG modeling is justified, we have performed QMC simulations on a 3D array of Heisenberg spin chains coupled weakly by almost Ising-like FM interactions (see Sec. (S5) of the SI for details of the model). It has been shown that this model reproduces the phase diagram of YbAlO 3 and exhibits three distinct ordered phases: low-field commensurate AFM, longitudinal SDW, and transverse canted AFM 21,26 , although it does not reproduce the M S /3 plateau (the plateau was not the primary focus of the present work, and we use DMRG to discuss its possible origin in Sec. (S6) of the SI). Figure 3d shows the calculated longitudinal spin structure factors, S zz (q, ω = 0), within the SDW phase at several representative magnetic fields. One can see that the signal consists of four peaks at q 1 = π δ, q 2 = 2π − 2δ, q 3 = π + 3δ, and q 4 = 2π−4δ, or in general nq 1 . Note that the intensity of the fundamental peak at q 1 significantly exceeds the weaker satellites at q 2 and q 3 , in good agreement with both the DMRG results and the experimental observations. Therefore, by combining the results of bosonization, DMRG and QMC we can conclude that the Umklapp scattering acts as an additional incommensurately modulated Zeeman potential of the form J ab cosð2k F rÞ in a single chain. This potential modifies the fermion bandstructure, as shown in Fig. 1a, and allows the multiple coherent scattering of fermions represented schematically in Fig. 3b, which gives rise to the observed satellites. In real space, these processes distort the SDW towards the more square-like shape shown in Fig. 3f1. It is worth noting that this mechanism for stabilizing higher-order harmonics should be more prominent with larger intrachain Ising anisotropy (Δ > 1), and thus one would expect that the satellites also appear in BaCo 2 V 2 O 8 and SrCo 2 V 2 O 8 15,16,27 . However, because S zz ðnq 1 Þ / ðJ ab =JÞ n , it may be difficult to resolve the satellites in ENS data because of the weak J ab /J ratio (≲ 1/30) in these compounds.
The physical consequences of multiple fermion scattering should manifest not only in reciprocal space-sensitive experiments, such as neutron diffraction, but also in low-temperature transport and thermodynamic properties. Umklapp scattering is known to suppress the thermal and electrical conductivity in metals, and thus multiple fermion scattering should provide an additional scattering channel reducing the thermal conductivity mediated by spin excitations. The presence of the incommensurate internal field in YbAlO 3 that allows for multiple scattering also opens a small gap in the spin excitation spectra, thus the spinon thermal conductivity of YbAlO 3 within the SDW phase would be altered in two ways. This interplay between the gapping effect and the extra contribution to the fermionic component of multiple scattering makes thermal transport measurements on YbAlO 3 highly desirable. The presence of the gap should also be observable by specific-heat measurements, which are similarly sensitive to the density of states of the magnetic quasiparticles. However, we note that a direct experimental detection in either quantity may be difficult, because the strength of multiple scattering scales with (J ab /J), and therefore is rather weak.
To summarize, we have observed distorted SDW ordering in YbAlO 3 by means of neutron scattering and have provided a detailed theoretical description. These results are of interest because we show how diffraction can be used to probe the bandstructure of fermionic quasiparticles in quantum magnets. We have demonstrated that the incommensurate effective field modulates the chemical potential and produces additional shadow bands as shown in Fig. 3b, which allow multiple (Umklapp) scattering of the fermions. To our knowledge, the coherent multiple scattering of fermions has not been detected before in either quantum magnets or one-dimensional metals. Thus, our results provide valuable insight into coherent quantum manybody phenomena and motivate further efforts to search for this effect in low-dimensional charge-and SDW condensed-matter systems, as well as in ultracold atoms in optical lattices.

Methods
Experimental details. A single crystal of YbAlO 3 with a mass of ≈ 0.5 g was prepared by a Czochralski technique, as described elsewhere 28,29 . The INS measurements were performed at the time-of-flight spectrometer CNCS 30,31 at the Spallation Neutron Source (SNS) at Oak Ridge National Laboratory with E i = 1.55 and 3.32 meV. To generate the magnetic field we used an 8 T vertical magnet. ENS measurements were performed using the TAS FLEXX at the Helmholtz Zentrum Berlin 25 with fixed k i = k f = 1.3 Å −1 . To increase the q-resolution we installed a 60' Söller collimator between the sample and the analyzer, and set both monochromator and analyzer with no horizontal focusing. The magnetic field was applied using the vertical cryomagnet VM-4. In both the neutron scattering experiments, the sample was oriented in the (0KL) scattering plane and the magnetic field was applied along with the easy a axis. The magnetization was measured using a high-resolution Faraday magnetometer 32 with a dilution cryostat. Reduction of the TOF data was performed using the Mantid 33 and Horace 34 software packages.
QMC and DMRG calculations. QMC simulations of the coupled spin-chain model were performed with a maximum system size of 20 × 20 × 128 sites and the lowest temperature accessed was T/J = 0.01. DMRG calculations were performed for two different spin models. The first one is given by Eq. (3) for the calculation of S zz (q, ω = 0) we used the standard DMRG method with open boundary conditions. We studied a 1D chain of L = 240 and 1600 density-matrix eigenstates that were kept in the renormalization procedure. In this way, the maximum truncation error, i.e., the discarded weight, can be negligible (≤10 −16 ). The second model was a two-leg spin-ladder with a next-neighbor frustrating intrachain interaction, which was used to provide a qualitative description of the M S /3 plateau, see Sec. (S6) of SI.