Intertwined density waves in a metallic nickelate

Nickelates are a rich class of materials, ranging from insulating magnets to superconductors. But for stoichiometric materials, insulating behavior is the norm, as for most late transition metal oxides. Notable exceptions are the 3D perovskite LaNiO3, an unconventional paramagnetic metal, and the layered Ruddlesden-Popper phases R4Ni3O10, (R = La, Pr, Nd). The latter are particularly intriguing because they exhibit an unusual metal-to-metal transition. Here, we demonstrate that this transition results from an incommensurate density wave with both charge and magnetic character that lies closer in its behavior to the metallic density wave seen in chromium metal than the insulating stripes typically found in single-layer nickelates like La2-xSrxNiO4. We identify these intertwined density waves as being Fermi surface-driven, revealing a novel ordering mechanism in this nickelate that reflects a coupling among charge, spin, and lattice degrees of freedom that differs not only from the single-layer materials, but from the 3D perovskites as well.


Crystal Structure of La4Ni3O10
The unit cell of La4Ni3O10 contains two trilayer perovskite slabs of corner-sharing NiO6 octahedra separated by a rocksalt layer of LaO (see Fig. 2a of main text). Depending on details of the growth conditions and post-growth treatment 1 , La4Ni3O10 crystallizes in either an orthorhombic (Bmab) or monoclinic (P21/a) crystal structure above TMMT. The crystal studied here shows TMMT = 140 K, indicating that the high temperature structure belongs to the monoclinic space group P21/a. However, this structure is pseudo-orthorhombic, and for ease of notation we employ the Bmab setting with lattice parameters of a~5.41 Å, b~5.46 Å, and c~27.97 Å throughout this report. In this setting, the two trilayers are symmetry-related by the B-centering operation.
Considering the limited pixels for each peak from linecuts from 15-ID-D data, we further performed single crystal diffuse scattering experiments at Beamline 33-BM-C at the Advanced Photon Source using a point detector at 70 K (λ=0.7749 Å). hkl scans were performed using steps of 0.0025 reciprocal lattice unit (r.l.u.) for k scans and 0.01 r.l.u. for l scans. Superlattice peaks were observed and fit to the Gaussian function (see Supplementary Fig. 2a) Fig. 2b shows x-ray scans collected on three unique single crystal La4Ni3O10 specimens to establish the consistency of the charge q-vector. This consistency tends to rule out incommensurability driven by defects, disorder, or nonstoichiometry. Further evidence is presented in the main text that argues against oxygen nonstoichiometry based on the amount needed to shift the propagation vector from the commensurate qc=2/3 expected for the nominal hole concentration to the observed incommensurate position. We estimate an O vacancy concentration ~4 times larger than measured would be required.
Supplementary Fig. 2b Comparison of charge qc-vector in three unique specimens. Sample A was measured at APS Beamline 15-ID-D, Sample B was measured at APS Beamline 33-BM-C, and Sample C was measured at APS Beamline 30-ID-C. The commensurate qc=1/3 position expected based on charge balance in a real-space stripe is marked.

Correlation length of CDW
Using data collected at the 33-BM-C beamline of the Advanced Photon Source, the correlation length of the charge order modulation was estimated through analysis of the peak width of the superlattice reflections. The correlation length x was determined as x =1/G, where G is the halfwidth at half-maximum of the Lorentzian component. A Gaussian fit was performed on the main Bragg peak (1,4,18) in order to estimate the instrumental resolution, and this Gaussian width was fixed in the Voigt function (a convolution of Gaussian and Lorentzian functions) used to fit the superlattice peaks ( Supplementary Fig. 3).
The in-plane and out-of-plane correlation length are xab=16-22b (87-118 Å) and xc=0.67-0.87c (18.8-24 Å), respectively. The estimated correlation length along c is comparable to the 21 Å that was found independently from 15-ID-D data shown in Supplementary Fig. 4. Supplementary Fig. 3 Correlation length of the CDW estimated from data collected at Beamline 33-BM-C at the Advanced Photon Source. data). The correlation length x was determined as x =1/G, where G is the half-width at halfmaximum of the Lorentzian component. A Gaussian fit was performed on the main Bragg peak (0,4,0) in order to estimate the instrumental resolution, and this Gaussian width was fixed in the Voigt function (a convolution of Gaussian and Lorentzian functions) used to fit the superlattice peak (-1,4.24,0). Using this procedure, G was determined to be 0.215 r.l.u. (where 1 r.l.u.=2p/c) yielding x= 0.7c=21 Å.

Neutron scattering in the hk0 plane
In the main text, neutron scattering measurements were shown in the 0kl plane. Here, we show measurements performed in the hk0 plane on the same instrument. Please note that whereas the measurements in the 0kl plane were performed with one specimen, the measurements in the hk0 plane were performed with four co-aligned specimens. A 180 K data set ( Supplementary Fig. 5b) has been subtracted from the 1.6 K data set ( Supplementary Fig. 5a) to reveal the magnetic superlattice peaks from the noise in the raw data. Supplementary Fig. 5c shows the map with the high temperature subtracted as a background. Supplementary Fig. 5d shows a line-cut through the temperature subtracted data along h. Note that the signal in this plane is much weaker than the maximal signal in the 0kl plane because the strong maxima occur at l=2 and 6, rather than 0. Nevertheless, a contribution to the hk0 plane is observed that may result from intrinsically broad peaks in l that are out of plane, exacerbated by the course out-of-plane resolution of the spectrometer.   Supplementary Fig. 7 Intensities for different assumed SDW stackings along the c-axis for the (0, 5/8, l) line cut, calculated using the equation listed in the main text. The titles indicate the stacking pattern along c (the cz weights of the main text). The upper left plot is the SDW model presented in the main text, which has prominent peaks at l=2 and 6.   0 2 4 6 8 10 12 14 16 (1,0,-1; 1,0,-1) Fig. 9 Intensities for different assumed CDW stackings along the c-axis for the (-1, 4.25, l) line cut. The titles represent the stacking patterns along c (numbers in parentheses are the cz weights described in the main text). The upper left plot is the CDW model presented in the main text, which has prominent peaks at l=0 and 14.  Above we have presented a simple model that accounts qualitatively for the distribution of intensity in the charge peaks measured by XRD. There are certain systematic discrepancies in the intensity distributions that highlight the need for a more sophisticated model. For example, as shown in Supplementary Fig. 11, our model for the (-1,k,0) cut of main text Fig. 2d finds strong peaks at k±qc (k odd) in agreement with the data. For the (0,k,19) cut of main text Fig. 2e, there are significant differences. In particular, our model indicates strong peaks (k even) and weak peaks (k odd), whereas the data find these peaks to be of comparable intensity. Moreover, the SL intensities are reduced relative to those in main text Fig. 2d by 10 2 as compared to about 10 1 found in experiment. To successfully reproduce the data, one would need a more sophisticated model than the present one. For instance, taking into account the monoclinic distortion is non-trivial in the model used here since b ≠ 90 degrees, meaning the simple cosine function would be de-phased when going from one layer to the next. Moreover, for hard x-rays, one is most sensitive to the atomic displacements. A strain wave model based on displacements of La, Ni, and O ions is a natural next step, albeit an involved endeavor, as can be seen from studies in the cuprates 2 . Supplementary Fig. 11 Intensities from the CDW model presented in the main text for the two k cuts shown in Fig. 2. For the right plot, note the alternation of strong peaks (k even) with weak peaks (k odd), in contrast to the data where these peaks have comparable intensity.

Real space spin-and charge-stripe model calculations
Various real space models for the density waves can be constructed by assuming 16 Ni ions per plane per magnetic unit cell. With an appropriate pattern of spins (d 7 ) and holes (d 8 ), the primary peak at 1-qs=5/8 can be generated, but because of the real space (square wave) form, (odd) harmonics appear with the strongest at 2-3qs=7/8, albeit with an intensity reduced by an order of magnitude relative to 1-qs =5/8 ( Supplementary Fig. 12, top). These harmonics are even larger (~17% of the fundamental) for the charge-stripe ( Supplementary Fig. 12, bottom left) and would not escape our observation in the data along (-1,k,0) shown in Fig. 2d of the main text. Given the absence of these harmonics in the data and the fact that the real space pattern does not correspond to the correct doping relative to LaNiO3 (which would imply a 1-qs=2/3 instead), we did not pursue real space models further. Supplementary Fig. 12 Intensities from a real space stripe model. For the top row, the title of the plots shows the SDW weights on each Ni ion for eight successive diagonal rows (±1 corresponding to d 7 , 0 to d 8 ). For the eight rows beyond this, the spins are reversed, with the pattern repeated. The stacking along c is the same as before (↑, −, ↓), with the lower and upper trilayers of a given unit cell being in-phase. This acts to maintain peaks at even l. Since there are five non-zero weights for each eight rows, this gives rise to a primary peak at 1-qs=5/8. The secondary peaks (k=7/8, etc.) are odd harmonics due to the "square wave" nature of the stripe, as compared to the sinusoidal behavior of the itinerant model. For the bottom row, the title of the plots shows the CDW weights on each Ni ion for 8 diagonal rows (d 7 being -1/3, d 8 being +2/3, relative to the nominal valence).
As with the itinerant model, the CDW is present on all planes, is in-phase between the three layers, and is out-of-phase between the trilayers. In both spin and charge cases, the BCT shift is taken into account by shifting the pattern by one place when going from the lower to the upper trilayer. Note for the charge stripe, the harmonics in the plane (k=3/2, etc.) are even larger than the spin harmonics (intensity is roughly 17% of the fundamentals at k=7/4, etc.).

Why the CDW is not on the oxygen sites
Our model assumes the density wave is centered on the Ni sites. This is sensible for the SDW, but what about for the CDW? A simple argument shows that the CDW cannot involve the planar oxygen sites. The in-plane planar oxygen coordinates (Bmab) are (¼,¼), (¾,¼), (¼,¾), and (¾,¾). For the first pair at y=1/4, the CDW will have the same amplitude, and similarly for the second pair at y=3/4. Therefore, for h odd, the sum over a pair vanishes because of the ! • structure factor, since their x coordinates differ by ½. This is in clear contradiction to experiment. If instead one assumes a d-wave form factor for the CDW, where the CDW amplitude has opposite phase on the two members of the pair, then the sum vanishes for h even, again inconsistent with the data. We argue that the small monoclinic distortion found in La4Ni3O10 1 can only cause weak violations of these rules. Finally, we remark that high energy x-rays are primarily sensitive to atomic displacements. Not only would such a strain wave differ from the CDW itself, but also other ions, such as the apical oxygens and La, can contribute, as seen in cuprates 2 . More information than we have at present would be needed to model such a complicated strain wave.

Band structure and Fermi surface nesting
Given the above context, we can ask whether Fermi surface nesting is relevant to the physics of La4Ni3O10, as demonstrated in chromium 3 and suggested as well for LaNiO3 4 . A recent ARPES study on La4Ni3O10 5 suggests that the opening of a gap in the 3z 2 -r 2 band accompanies the CDW at TMMT. We note that La4Ni3O10 crystallizes in an orthorhombic (Bmab) or monoclinic (P21/a) supercell of the tetragonal I4/mmm parent, leading to backfolding of the bands. This is evident in Supplementary Fig. 13a, where the backfolding occurs exactly where Ref. 5 asserts the CDW gap to lie. We conjecture that this backfolding is likely the gap that the ARPES is observing. Higher resolution ARPES data will likely be needed to detect the actual CDW gap.
We note that Seo et al. 6 have proposed from extended Hückel tight-binding calculations that the MMT arises from a CDW instability resulting from Fermi surface nesting. To explore this possibility further, we performed DFT calculations for La4Ni3O10 using Bmab coordinates. We then fit the eight bands nearest to the Fermi energy with a Fourier series spline fit. The fit is shown in Supplementary Fig. 13a and the resulting Fermi surfaces in the kz=0 plane in Supplementary  Fig. 13b. Note the zone folding relative to a hypothetical undistorted I4/mmm structure. Since the bands 1 and 2 Fermi surfaces are hole-like, and the bands 3 and 4 Fermi surfaces are electron-like, nesting is a distinct possibility. To test this possibility, we first accurately determined the Fermi energy and then calculated the static susceptibility (c) using a linear tetrahedron method. The resulting interband c terms are shown in Supplementary Fig. 13c. Here we present a cut along (q,0,0). Almost equivalent results are found along (0,q,0), implying that the orientation of the density wave with ∥ is due to pinning by the structure. The band 2 to band 3 interband term shows a maximum at (0.5,0,0) but has subsidiary peaks at (0.5±0.11,0,0), with 0.61 lying close to the observed 1-qs=0.62.
Since we find considerable structure in the density of states near the Fermi energy, we recalculated c by reducing the Fermi energy by 30 meV (Supplementary Fig. 13d). Now the side peaks become the global maxima, although they are shifted to (0.5±0.09,0,0). These results indicate that nesting is a distinct possibility. To improve on these results would require more sophisticated calculations that include not only spin matrix elements between the various band eigenvectors, but also potential correlation-induced energy shifts between the Ni 3d and O 2p states. For these reasons, we did not pursue further calculations along these lines, including the qz dependence (which we find to be weak). Here, q is in 2p/a units. The outermost of the three peaks in c for the 2-3 interband term is at q=0.61, close to the observed SDW wavevector at 1-qs=0.62. This peak becomes a global maximum in d, though it is shifted to q=0.59. This indicates that nesting is a potential cause of the SDW in La4Ni3O10.