^†^†thanks: This manuscript has been authored by UT-Battelle, LLC under Contract No. DE-AC05-00OR22725 with the U.S. Department of Energy. The United States Government retains and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a non-exclusive, paid-up, irrevocable, world-wide license to publish or reproduce the published form of this manuscript, or allow others to do so, for United States Government purposes. The Department of Energy will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan (http://energy.gov/downloads/doe-public-access-plan).

Codimension-Two Spiral Spin-Liquid in the Effective Honeycomb-Lattice Compound Cs₃Fe₂Cl₉

Shang Gao sgao@ustc.edu.cn Department of Physics, University of Science and Technology of China, Hefei, Anhui 230026, People’s Republic of China Materials Science & Technology Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Chris Pasco Materials Science & Technology Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Otkur Omar Department of Physics, University of Science and Technology of China, Hefei, Anhui 230026, People’s Republic of China Qiang Zhang Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Daniel M. Pajerowski Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Feng Ye Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Matthias Frontzek Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Andrew F. May Materials Science & Technology Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Matthew B. Stone stonemb@ornl.gov Neutron Scattering Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA Andrew D. Christianson christiansad@ornl.gov Materials Science & Technology Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA

(May 29, 2024)

Abstract

A codimension-two spiral spin-liquid is a correlated paramagnetic state with one-dimensional ground state degeneracy hosted within a three-dimensional lattice. Here, via neutron scattering experiments and numerical simulations, we establish the existence of a codimension-two spiral spin-liquid in the effective honeycomb-lattice compound Cs₃Fe₂Cl₉ and demonstrate the selective visibility of the spiral surface through phase tuning. In the long-range ordered regime, competing spiral and spin density wave orders emerge as a function of applied magnetic field, among which a possible order-by-disorder transition is identified.

pacs:

Introduction. A spiral spin-liquid (SSL) is an exotic type of correlated paramagnetic state where the low energy dynamics consist of collective spiral correlations Bergman et al. (2007); Lee and Balents (2008); Mulder et al. (2010); Zhang and Lamas (2013); Niggemann et al. (2020); Attig and Trebst (2017); Balla et al. (2020, 2019); Lee and Balents (2008); Chen (2017); Yao et al. (2021); Huang et al. (2022); Buessen et al. (2018); Liu et al. (2020); Mohylna et al. (2022). A characteristic feature of a SSL is that the propagation vectors of the degenerate spiral ground states form a continuous surface in reciprocal space Bergman et al. (2007). Such an unusual yet clearly defined feature has stimulated a strong interest by the community to experimentally identify and understand SSLs in real materials Tristan et al. (2005); T. Suzuki and H. Nagai and M. Nohara and H. Takagi (2007); Krimmel et al. (2009); Matsuda et al. (2010); MacDougall et al. (2011); Zaharko et al. (2011); Nair et al. (2014); MacDougall et al. (2016); Ge et al. (2017); Graham et al. (2023); Gao et al. (2017, 2020); Rosales et al. (2022); Guratinder et al. (2022); Chamorro et al. (2018); Bai et al. (2019); Tsurkan et al. (2021); Haraguchi et al. (2019); Abdeldaim et al. (2020); Otsuka et al. (2020); Wessler et al. (2020); Bordelon et al. (2021); Balz et al. (2016); Pohle et al. (2021); Takahashi et al. (2024); Gao et al. (2022a); Cole et al. (2023); Gao et al. (2022b); Hsieh et al. (2022). Through experimental and theoretical studies, compounds with a bipartite lattice, e.g. the honeycomb Balz et al. (2016); Pohle et al. (2021); Gao et al. (2022a); Cole et al. (2023) and diamond Graham et al. (2023); Balz et al. (2016); Pohle et al. (2021); Gao et al. (2017); Bai et al. (2019); Guratinder et al. (2022) lattices, have been demonstrated as the most fertile hosts of SSLs.

Refer to caption — Figure 1: (a) AB-stacked triangular bilayers formed by the Fe³⁺ ions in Cs₃Fe₂Cl₉. Atoms belonging to the neighboring bilayers are shown in red and blue, respectively. The $J_{1}$ , $J_{2}$ , and $J_{3}$ bonds are shown by thick black lines, thin black lines, and thin gray lines, respectively. Arrows indicate the spin directions of the collinear ground state with $\bm{q}=(\frac{1}{2},0,0)$ . (b) The AB-stacked triangular bilayers viewed along the $c$ axis. (c) Refinement result of the powder neutron diffraction data measured on POWGEN at $T$ = 1.6 K. Data points are shown as red circles. The calculated pattern is shown as the black solid line. The vertical bars indicate the positions of the structural (upper) and magnetic (lower) Bragg peaks for Cs₃Fe₂Cl₉. The blue line at the bottom shows the difference of measured and calculated intensities. The goodness-of-fit parameters are $R_{\textrm{p}}=18.6\%$ and $R_{\textrm{wp}}=10.2\%$ .

Identifying novel SSL hosts is crucial for the realization of exotic spin textures like skyrmions Gao et al. (2020); Rosales et al. (2022); Takeda et al. (2024) and subdimensional quasiparticles like fractons Yan and Reuther (2022); Pretko et al. (2020), and will also establish new candidate compounds to study the thermal and quantum order-by-disorder (ObD) transitions that are elusive in real materials Villain et al. (1980); Henley (1989); Green et al. (2018); Bergman et al. (2007); Mulder et al. (2010). According to theoretical studies, SSLs can be classified by their codimension, a quantity that characterizes the dimensional difference between the spiral surface and the host system Yao et al. (2021). Experimentally identified SSLs, including those observed on the diamond and honeycomb lattices, exhibit either a two-dimensional (2D) spiral surface on a three-dimensional (3D) lattice Gao et al. (2017); Bai et al. (2019); Guratinder et al. (2022); Gao et al. (2022b); Graham et al. (2023) or a one-dimensional (1D) spiral surface, i.e. a degenerate line, on a 2D lattice Balz et al. (2016); Pohle et al. (2021); Takahashi et al. (2024); Gao et al. (2022a), thus all falling in the codimension one category. Although codimension-two SSLs have been predicted to exist on the AB-stacked triangular lattice Yao et al. (2021); Liu et al. (2016); Hoang and Diep (2012); Rastelli et al. (1983), their experimental realizations remain an open question.

Recent transport and magnetic characterizations of Cs₃Fe₂Cl₉ unveil possible honeycomb physics on a 3D lattice Ishii et al. (2021). In this compound, magnetic Fe³⁺ ions with spin $S=\frac{5}{2}$ form AB-stacked triangular bilayers as shown in Figs. 1(a) and (b). Under dominant ferromagnetic (FM) interactions with exchange energy strength $J_{1}$ within the bilayers, a codimension-two SSL with a 1D spiral surface in the integer- $l$ planes may exist for $|J_{3}/J_{2}|>1/6$ , where the threshold ratio is the same as that of the codimension-one SSL on the original honeycomb lattice Mulder et al. (2010). Although the magnetic ground state and spin dynamics in Cs₃Fe₂Cl₉ remain unexplored, a rich phase diagram has been established Ishii et al. (2021), suggesting the existence of strong magnetic frustration in this compound.

Here, through elastic and inelastic neutron scattering experiments on both single crystal and polycrystalline samples, we show that a codimension-two SSL with a uniaxial anisotropy is realized in Cs₃Fe₂Cl₉. The sprial surface’s visibility is phase tuned as a function of the wavevector normal to the honeycomb plane. By combining neutron diffraction experiments and classical Monte Carlo simulations, we clarify the eight field-induced ordered phases as competing spiral and spin density wave (SDW) orders, among which a possible order-by-disorder transition is identified.

Magnetic ground state. Powder neutron diffraction experiments were performed on POWGEN Huq et al. (2011) at the Spallation Neutron Source (SNS) of the Oak Ridge National Laboratory (ORNL) to determine the magnetic LRO in Cs₃Fe₂Cl₉ below $T_{N}\sim 5.4$ K. Details on the sample preparation and neutron scattering experiments are presented in the Supplemental Material sup . As shown in Fig. 1(c), magnetic Bragg peaks belonging to the propagation vector $\bm{q}_{\textrm{I}}=(\frac{1}{2},0,0)$ are observed at low temperatures. Through Rietveld refinements Rodriguez-Carvajal (1993), the magnetic ground state is determined to be collinear as shown by arrows in Fig. 1(a), with an ordered moment size of 4.23(6) $\mu_{\textrm{B}}$ . This magnetic order is similar to the magnetic ground state of the isostructural Cs₃Fe₂Br₉ Brüning et al. (2021).

Spin Dynamics and Modeling. To determine the exchange coupling strengths, inelastic neutron scattering (INS) experiments were performed at CNCS Ehlers et al. (2011) at SNS of ORNL sup . Figure 2(a) presents the INS spectra collected at $T=2$ K with an incident neutron energy of $E_{i}=3.32$ meV. Even for the powder sample, it is clear that there are two highly dispersive magnon modes centered around $E\sim 0.6$ and 1.8 meV energy transfer. From Fig. 2(a) and also the $E_{i}=1.0$ meV data in the Supplemental Material sup , an excitation gap of $\Delta\sim 0.2$ meV is observed at wavevector transfer $Q\sim 0.6$ Å^-1, suggesting the existence of a uniaxial single-ion anisotropy (SIA) that stabilizes the collinear ground state.

Using linear spin wave theory as implemented in the SpinW program Toth and Lake (2015), $\chi^{2}$ -fits to the INS spectra were performed to analyze the spin interactions. A minimal $J_{123}$ - $D_{z}$ model with Hamiltonian $\mathcal{H}=\sum_{\langle ij\rangle\in n}J_{n}\bm{S}_{i}\cdot\bm{S}_{j}+D_{z}(% S_{z})^{2}$ is employed in our calculations, which considers Heisenberg exchange interactions up to the third neighbors as shown in Fig. 1 plus a uniaxial SIA term, $D_{z}$ . As compared in Fig. 2, the $J_{123}$ - $D_{z}$ model with fitted coupling strengths of $J_{1}=-0.238(2)$ , $J_{2}=0.085(2)$ , $J_{3}=0.053(1)$ , and $D_{z}=-0.037(2)$ meV reproduces the INS spectra. The dominant FM $J_{1}$ together with comparable strengths of $J_{2}$ and $J_{3}$ favor SSLs Mulder et al. (2010); Yao et al. (2021); Liu et al. (2016), while the relatively high magnitude of $|D_{z}|\sim 0.65J_{3}$ indicates the importance of the SIA.

Codimension-two SSLs. Figure 3(a) presents the diffuse neutron scattering pattern for Cs₃Fe₂Cl₉ in the ( $h$ , $k$ , 0) plane. Data were collected from CORELLI Ye et al. (2018) at SNS of ORNL using an $\sim 8$ mg crystal sup . At $T=6$ K, triangular shaped lobes are observed around the $K$ -( $\frac{1}{3}$ , $\frac{1}{3}$ , 0) points. As indicated by the dashed curves, the shape of the spiral surface is reproduced by a $J_{1}^{h}$ - $J_{2}^{h}$ honeycomb-lattice model with $J_{2}^{h}/J_{1}^{h}=0.62$ , where $J_{1}^{h}$ $(J_{2}^{h})$ is equivalent to $J_{2}$ ( $J_{3}$ ) in the $J_{123}$ - $D_{z}$ model. The equal-time spin correlations for the fitted $J_{123}$ - $D_{z}$ model can be calculated using the self-consistent Gaussian approximation (SCGA) method Canals and Garanin (2001); sup . Since the critical correlations are underestimated in the SCGA method, a reduced $T=5$ K is assumed in the calculations to better describe the experimental data sup . As compared in Fig. 3(a), the calculated and experimental intensity distributions agree well with each other, both following the spiral surface of a SSL on a honeycomb lattice.

A special feature of the SSL in Cs₃Fe₂Cl₉ is the visibility of the spiral surface, which arises from its unique codimension two. On the original honeycomb lattice with antiferromagnetic $J_{1}^{h}$ , the spiral surface within the first Brillouin zone has a structure factor of zero due to the interference between the two sublattices Liu et al. (2016), leading to a diffuse pattern similar that shown in Fig. 3(a). However, as compared in Figs. 3(a) and (b), the visibility of the spiral surface in Cs₃Fe₂Cl₉ is complementary between the $l=0$ and $l=1$ planes. This variance in visibility originates from the $l$ -dependence of the phase factor, which modulates the interference between the two sublattices and thus reflects the higher codimension of the SSL Liu et al. (2016). By integrating the intensity in the range of $l=[-1.1,1.1]$ r.l.u., the complete spiral surface can be recovered in Fig. 3(c). The scattering pattern in the ( $h,~{}h,~{}l$ ) plane presented in Fig. 3(d) reveals weak scattering intensity between the integer- $l$ planes. According to our calculations, the intensity of the interplanar scattering diminishes with decreasing temperatures, thus can be attributed to thermal excitations out of the ground state manifolds. This observation is further confirmed in Fig. 3(e) through the measured and calculated $l$ -dependence of the integrated intensity within a rectangular area of $[0.125,0.5]$ r.l.u. and $[-0.125,0.125]$ r.l.u. along the $(h,h,0)$ and $(h,-h,0)$ directions, respectively. In the same area, the center of mass for the scattering intensity, shown in Fig. 3(f), varies continuously along the $(h,h,0)$ direction in excellent agreement with our model.

Competing orders and a possible ObD transition. In the absence of SIA, SSLs on a honeycomb lattice have been predicted to exhibit an ObD transition at low temperatures Bergman et al. (2007); Mulder et al. (2010); Liu et al. (2020). In the regime of $J_{2}^{h}/J_{1}^{h}>0.5$ , the magnetic propagation vector $\bm{q}$ of the ObD phase lies at the corners of the triangular-shaped spiral surface Mulder et al. (2010); Liu et al. (2020). This is obviously not the case for Cs₃Fe₂Cl₉ with $\bm{q}_{\textrm{I}}=(\frac{1}{2},0,0)$ as revealed in Fig. 1(c), suggesting the strong impact of the SIA.

However, the rich phase diagram reported for Cs₃Fe₂Cl₉ Ishii et al. (2021), which is also reproduced in Fig. 4(a), indicates the the possibility of an ObD in magnetic field. In fields of $3\lesssim H\lesssim 6$ T, multiple phases, II-VI, emerge. Two additional phases exist at $H\gtrsim 6$ T, including a $\frac{1}{2}$ -magnetization plateau phase, VII, up to $\sim 12$ T and a transitional phase, VIII, that precedes the field-polarized FM phase at $H\gtrsim 15$ T.

Single crystal neutron diffraction experiments were performed on WAND² Frontzek et al. (2018) at the High Flux Isotope Reactor (HFIR) at ORNL to clarify the LRO phases as a function of magnetic field up to 6 T sup . As summarized in Figs. 4(b)-4(e), four different types of diffraction patterns are observed in the ( $h$ , $k$ , 0) plane. Compared to that in phase I, the diffraction pattern of phase II in Fig. 4(c) exhibits additional weak magnetic Bragg peaks at $\bm{q}_{\textrm{II}}=(\frac{1}{4},0,0)$ , indicating the coexistence of minor $\uparrow\uparrow\uparrow\downarrow$ and major $\uparrow\downarrow\uparrow\downarrow$ magnetic domains, where $\uparrow$ ( $\downarrow$ ) represents spins (antiparallel) parallel with $c$ . Surprisingly, the diffraction patterns in phases III-V are similar to each other sup : all exhibiting main reflections at $\bm{q}_{\textrm{III}}=\frac{5}{6}\times(\frac{1}{3},\frac{1}{3},0)$ close to the center of the edges of the triangular-shaped spiral surface as indicated by dashed lines, together with weaker secondary reflections at $2\bm{q}_{\textrm{III}}$ . The existence of the latter often suggests a field-distorted spiral or SDW phase Gignoux et al. (1998); Stüßer et al. (2002). In phase VI, magnetic reflections belonging to $\bm{q}_{\textrm{VI}}=\frac{11}{9}\times(\frac{1}{3},\frac{1}{3},0)$ appear at the corners of the triangular-shaped spiral surface, which is exactly the position predicted by the ObD transition Mulder et al. (2010); Liu et al. (2020).

Due to the limited number of magnetic reflections collected on WAND², we resort to classical Monte Carlo simulations to understand the rich phases in Cs₃Fe₂Cl₉. As explained in the Supplemental Material sup , the fitted $J_{123}$ - $D_{z}$ model exhibits a relatively simple phase diagram that contains only two main phases with $\bm{q}_{\textrm{I}}=(\frac{1}{2},0,0)$ . Therefore, starting from the fitted $J_{123}$ - $D_{z}$ model, we examined the impacts of weak perturbations from isotropic exchange interactions up to the tenth neighbors one-by-one. The exchange paths of each perturbations and the representative theoretical phase diagrams are presented in the Supplemental Material sup .

Following this method, we found that by adding a weak fifth-neighbor coupling $J_{5}=-0.008$ meV that only marginally impacts the INS spectra and diffuse scattering patterns sup , the calculated phase diagram, shown in Fig. 4(f), reproduces the main phases observed in experiments. Phases VI^′ and VI^′′ (Phase V^′) adjacent to phases I and VII exhibit similar diffraction pattern as that in phase VI (V). These primed phases, V^′, VI^′, and VI^′′, are absent in experiments and may arise from finite size effects in our simulations. The magnetic orders in the main phases are described in Figs. 4(g)-(j), with their corresponding Ansatz and structure factor plots presented in the Supplemental Material sup . As is consistent with the experimental data, similar propagation vectors are observed in phases III-V. Among them, phases III and IV are determined to be of the spiral-type orders by inspecting the spin configuration snapshots in the Monte Carlo simulations, where the transition between phases III and IV can be attributed to the modulation of the $\bm{q}=0$ magnitude sup . As a contrast, both phases V and VI are revealed to be of the SDW orders with spins along the $c$ axis despite their different $\bm{q}$ vectors. In phases III-V where sufficient numbers of magnetic Bragg peaks can be extracted in our neutron diffraction experiments, the proposed magnetic structures agree with the experimental dataset sup , thus verifying the results of our Monte Carlo simulations.

Since $\bm{q}_{\textrm{VI}}$ in phase VI is the predicted propagation vector of the ObD theory Mulder et al. (2010); Liu et al. (2020), it is tantalizing to ascribe the transition from the SSL state to phase VI as entropy-driven. Compared to the theoretical phase diagram for the fitted $J_{123}$ - $D_{z}$ model, it can be concluded that in the presence of a magnetic field, the $J_{5}$ perturbations favor orders with propagation vectors over the spiral surface, among which $\bm{q}_{\textrm{VI}}$ is stabilized by thermal fluctuations at higher temperatures. The proposed SDW character of phase VI, however, poses a challenge in accurately determining the free energy for states above $T_{N}$ under the current theoretical framework Bergman et al. (2007), making it difficult to exclude perturbations other than thermal fluctuations.

Conclusion. Our neutron scattering experiments establish Cs₃Fe₂Cl₉ as a host of codimension-two SSLs, where the high codimension is manifested through the phase tuning of the visibility of the spiral surface. By combining neutron diffraction and classical Monte Carlo simulations, we clarify the rich field-induced phases as competing spiral and SDW orders, revealing Cs₃Fe₂Cl₉ as a candidate compound to study the ObD transition in frustrated magnets.

Acknowledgements.

We acknowledge helpful discussions with Gang Chen, Xu-Ping Yao, and James Jun He. This work was supported by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division. This research used resources at the Spallation Neutron Source (SNS) and the High Flux Isotope Reactor (HFIR), both are DOE Office of Science User Facilities operated by the Oak Ridge National Laboratory (ORNL). The proposal numbers are IPTS-28988 (POWGEN), 29069 (CNCS), 31251 (CORELLI), and 29937 (WAND²). Numerical calculations in this paper were partly completed on the supercomputing system in the Supercomputing Center of USTC. Works at USTC were funded by the National Science Foundation of China (NSFC) under the Grant No. 12374152.

References

Bergman et al. (2007) D. Bergman, J. Alicea, E. Gull, S. Trebst, and L. Balents, Nat. Phys. 3, 487 (2007).
Lee and Balents (2008) S. B. Lee and L. Balents, Phys. Rev. B 78, 144417 (2008).
Mulder et al. (2010) A. Mulder, R. Ganesh, L. Capriotti, and A. Paramekanti, Phys. Rev. B 81, 214419 (2010).
Zhang and Lamas (2013) H. Zhang and C. A. Lamas, Phys. Rev. B 87, 024415 (2013).
Niggemann et al. (2020) N. Niggemann, M. Hering, and J. Reuther, J. Phys.: Condens. Matter 32, 024001 (2020).
Attig and Trebst (2017) J. Attig and S. Trebst, Phys. Rev. B 96, 085145 (2017).
Balla et al. (2020) P. Balla, Y. Iqbal, and K. Penc, Phys. Rev. Research 2, 043278 (2020).
Balla et al. (2019) P. Balla, Y. Iqbal, and K. Penc, Phys. Rev. B 100, 140402(R) (2019).
Chen (2017) G. Chen, Phys. Rev. B 96, 020412(R) (2017).
Yao et al. (2021) X.-P. Yao, J. Q. Liu, C.-J. Huang, X. Wang, and G. Chen, Frontiers Phys. 16, 53303 (2021).
Huang et al. (2022) C.-J. Huang, J. Q. Liu, and G. Chen, Phys. Rev. Research 4, 013121 (2022).
Buessen et al. (2018) F. L. Buessen, M. Hering, J. Reuther, and S. Trebst, Phys. Rev. Lett. 120, 057201 (2018).
Liu et al. (2020) J. Q. Liu, F.-Y. Li, G. Chen, and Z. Wang, Phys. Rev. Research 2, 033260 (2020).
Mohylna et al. (2022) M. Mohylna, F. A. Gómez Albarracín, M. Žukovič, and H. D. Rosales, Phys. Rev. B 106, 224406 (2022).
Tristan et al. (2005) N. Tristan, J. Hemberger, A. Krimmel, H.-A. Krug von Nidda, V. Tsurkan, and A. Loidl, Phys. Rev. B 72, 174404 (2005).
T. Suzuki and H. Nagai and M. Nohara and H. Takagi (2007) T. Suzuki and H. Nagai and M. Nohara and H. Takagi, J. Phys.: Condens. Matter 19, 145265 (2007).
Krimmel et al. (2009) A. Krimmel, H. Mutka, M. M. Koza, V. Tsurkan, and A. Loidl, Phys. Rev. B 79, 134406 (2009).
Matsuda et al. (2010) M. Matsuda, M. Azuma, M. Tokunaga, Y. Shimakawa, and N. Kumada, Phys. Rev. Lett. 105, 187201 (2010).
MacDougall et al. (2011) G. J. MacDougall, D. Gout, J. L. Zarestky, G. Ehlers, A. Podlesnyak, M. A. McGuire, D. Mandrus, and S. E. Nagler, Proc. Natl. Acad. Sci. U.S.A. 108, 15693 (2011).
Zaharko et al. (2011) O. Zaharko, N. B. Christensen, A. Cervellino, V. Tsurkan, A. Maljuk, U. Stuhr, C. Niedermayer, F. Yokaichiya, D. N. Argyriou, M. Boehm, and A. Loidl, Phys. Rev. B 84, 094403 (2011).
Nair et al. (2014) H. S. Nair, Z. Fu, J. Voigt, Y. X. Su, and T. Brückel, Phys. Rev. B 89, 174431 (2014).
MacDougall et al. (2016) G. J. MacDougall, A. A. Aczel, Y. Su, W. Schweika, E. Faulhaber, A. Schneidewind, A. D. Christianson, J. L. Zarestky, H. D. Zhou, D. Mandrus, and S. E. Nagler, Phys. Rev. B 94, 184422 (2016).
Ge et al. (2017) L. Ge, J. Flynn, J. A. M. Paddison, M. B. Stone, S. Calder, M. A. Subramanian, A. P. Ramirez, and M. Mourigal, Phys. Rev. B 96, 064413 (2017).
Graham et al. (2023) J. N. Graham, N. Qureshi, C. Ritter, P. Manuel, A. R. Wildes, and L. Clark, Phys. Rev. Lett. 130, 166703 (2023).
Gao et al. (2017) S. Gao, O. Zaharko, V. Tsurkan, Y. Su, J. S. White, G. S. Tucker, B. Roessli, F. Bourdarot, R. Sibille, D. Chernyshov, T. Fennell, A. Loidl, and C. Rüegg, Nat. Phys. 13, 157 (2017).
Gao et al. (2020) S. Gao, H. D. Rosales, F. A. Gómez Albarracín, V. Tsurkan, G. Kaur, T. Fennell, P. Steffens, M. Boehm, P. Čermák, A. Schneidewind, E. Ressouche, D. C. Cabra, C. Rüegg, and O. Zaharko, Nature 586, 37 (2020).
Rosales et al. (2022) H. D. Rosales, F. A. G. Albarracín, K. Guratinder, V. Tsurkan, L. Prodan, E. Ressouche, and O. Zaharko, Phys. Rev. B 105, 224402 (2022).
Guratinder et al. (2022) K. Guratinder, V. Tsurkan, L. Prodan, L. Keller, J. P. Embs, F. Juranyi, M. Medarde, C. Rüegg, and O. Zaharko, Phys. Rev. B 105, 174423 (2022).
Chamorro et al. (2018) J. R. Chamorro, L. Ge, J. Flynn, M. A. Subramanian, M. Mourigal, and T. M. McQueen, Phys. Rev. Mater. 2, 034404 (2018).
Bai et al. (2019) X. Bai, J. A. M. Paddison, E. Kapit, S. M. Koohpayeh, J.-J. Wen, S. E. Dutton, A. T. Savici, A. I. Kolesnikov, G. E. Granroth, C. L. Broholm, J. T. Chalker, and M. Mourigal, Phys. Rev. Lett. 122, 097201 (2019).
Tsurkan et al. (2021) V. Tsurkan, H.-A. Krug von Nidda, J. Deisenhofer, P. Lunkenheimer, and A. Loidl, Phys. Rep. 926, 1 (2021).
Haraguchi et al. (2019) Y. Haraguchi, K. Nawa, C. Michioka, H. Ueda, A. Matsuo, K. Kindo, M. Avdeev, T. J. Sato, and K. Yoshimura, Phys. Rev. Mater. 3, 124406 (2019).
Abdeldaim et al. (2020) A. H. Abdeldaim, T. Li, L. Farrar, A. A. Tsirlin, W. Yao, A. S. Gibbs, P. Manuel, P. Lightfoot, G. J. Nilsen, and L. Clark, Phys. Rev. Mater. 4, 104414 (2020).
Otsuka et al. (2020) A. Otsuka, Y. Shimizu, G. Saito, M. Maesato, A. Kiswandhi, T. Hiramatsu, Y. Yoshida, H. Yamochi, M. Tsuchiizu, Y. Nakamura, H. Kishida, and H. Ito, Bull. Chem. Soc. Jpn. 93, 260 (2020).
Wessler et al. (2020) C. Wessler, B. Roessli, K. W. Krämer, B. Delley, O. Waldmann, L. Keller, D. Cheptiakov, H. B. Braun, and M. Kenzelmann, Npj Quantum Mater. 5, 85 (2020).
Bordelon et al. (2021) M. M. Bordelon, C. Liu, L. Posthuma, E. Kenney, M. J. Graf, N. P. Butch, A. Banerjee, S. Calder, L. Balents, and S. D. Wilson, Phys. Rev. B 103, 014420 (2021).
Balz et al. (2016) C. Balz, B. Lake, J. Reuther, H. Luetkens, R. Schonemann, T. Herrmannsdorfer, Y. Singh, A. T. M. Nazmul Islam, E. M. Wheeler, J. A. Rodriguez-Rivera, T. Guidi, G. G. Simeoni, C. Baines, and H. Ryll, Nat. Phys. 12, 942 (2016).
Pohle et al. (2021) R. Pohle, H. Yan, and N. Shannon, Phys. Rev. B 104, 024426 (2021).
Takahashi et al. (2024) H. Takahashi, C.-C. Hsu, F. Jerzembeck, J. Murphy, J. Ward, J. D. Enright, J. Knapp, P. Puphal, M. Isobe, Y. Matsumoto, H. Takagi, J. C. Séamus Davis, and S. J. Blundell, arXiv e-prints , arXiv:2405.02075 (2024), arXiv:2405.02075 [cond-mat.str-el] .
Gao et al. (2022a) S. Gao, M. A. McGuire, Y. Liu, D. L. Abernathy, C. d. Cruz, M. Frontzek, M. B. Stone, and A. D. Christianson, Phys. Rev. Lett. 128, 227201 (2022a).
Cole et al. (2023) A. Cole, A. Streeter, A. O. Fumega, X. Yao, Z.-C. Wang, E. Feng, H. Cao, J. L. Lado, S. E. Nagler, and F. Tafti, Phys. Rev. Mater. 7, 064401 (2023).
Gao et al. (2022b) S. Gao, G. Pokharel, A. F. May, J. A. M. Paddison, C. Pasco, Y. Liu, K. M. Taddei, S. Calder, D. G. Mandrus, M. B. Stone, and A. D. Christianson, Phys. Rev. Lett. 129, 237202 (2022b).
Hsieh et al. (2022) T.-C. Hsieh, H. Ma, and L. Radzihovsky, Phys. Rev. A 106, 023321 (2022).
Takeda et al. (2024) H. Takeda, M. Kawano, K. Tamura, M. Akazawa, J. Yan, T. Waki, H. Nakamura, K. Sato, Y. Narumi, M. Hagiwara, M. Yamashita, and C. Hotta, Nat. Commun. 15, 566 (2024).
Yan and Reuther (2022) H. Yan and J. Reuther, Phys. Rev. Research 4, 023175 (2022).
Pretko et al. (2020) M. Pretko, X. Chen, and Y. You, Int. J. Mod. Phys. A 35, 2030003 (2020).
Villain et al. (1980) J. Villain, R. Bidaux, J.-P. Carton, and R. Conte, Journal de Physique 41, 1263 (1980).
Henley (1989) C. L. Henley, Phys. Rev. Lett. 62, 2056 (1989).
Green et al. (2018) A. G. Green, G. Conduit, and F. Krüger, Ann. Rev. Condens. Matter Phys. 9, 59 (2018).
Liu et al. (2016) D. T. Liu, F. J. Burnell, L. D. C. Jaubert, and J. T. Chalker, Phys. Rev. B 94, 224413 (2016).
Hoang and Diep (2012) D.-T. Hoang and H. T. Diep, Phys. Rev. E 85, 041107 (2012).
Rastelli et al. (1983) E. Rastelli, L. Reatto, and A. Tassi, J. Phys. C: Solid State Phys. 16, L331 (1983).
Ishii et al. (2021) Y. Ishii, Y. Narumi, Y. Matsushita, M. Oda, T. Kida, M. Hagiwara, and H. K. Yoshida, Phys. Rev. B 103, 104433 (2021).
Huq et al. (2011) A. Huq, J. P. Hodges, O. Gourdon, and L. Heroux, Z. Kristallogr. 1, 127 (2011).
(55) See Supplemental Material for details on sample preparations and characterizations, neutron scattering experiments and data analysis, and classical Monte Carlo simulations.
Rodriguez-Carvajal (1993) J. Rodriguez-Carvajal, Physica B: Conden. Matter 192, 55 (1993).
Brüning et al. (2021) D. Brüning, T. Fröhlich, D. Gorkov, I. Císařová, Y. Skourski, L. Rossi, B. Bryant, S. Wiedmann, M. Meven, A. Ushakov, S. V. Streltsov, D. Khomskii, P. Becker, L. Bohatý, M. Braden, and T. Lorenz, Phys. Rev. B 104, 064418 (2021).
Ehlers et al. (2011) G. Ehlers, A. A. Podlesnyak, J. L. Niedziela, E. B. Iverson, and P. E. Sokol, Rev. Sci. Instrum. 82, 085108 (2011).
Toth and Lake (2015) S. Toth and B. Lake, J. Phys.: Condens. Matter 27, 166002 (2015).
Ye et al. (2018) F. Ye, Y. Liu, R. Whitfield, R. Osborn, and S. Rosenkranz, J. Appl. Cryst. 51, 315 (2018).
Canals and Garanin (2001) B. Canals and D. A. Garanin, Can. J. Phys. 79, 1323 (2001).
Frontzek et al. (2018) M. D. Frontzek, R. Whitfield, K. M. Andrews, A. B. Jones, M. Bobrek, K. Vodopivec, B. C. Chakoumakos, and J. A. Fernandez-Baca, Rev. Sci. Instrum. 89, 092801 (2018).
Gignoux et al. (1998) D. Gignoux, D. Schmitt, and T. Shigeoka, Europhys. Lett. 44, 374 (1998).
Stüßer et al. (2002) N. Stüßer, U. Schotte, A. Hoser, M. Meschke, M. Meißner, and J. Wosnitza, J. Phys: Condens. Matter 14, 5161 (2002).

Supplemental Materials for:
Codimension-Two Spiral Spin-Liquid in the Effective Honeycomb-Lattice Compound Cs₃Fe₂Cl₉

I Powder sample synthesis

Cs₃Fe₂Cl₉ is extremely air sensitive and thus all handling of reagents and products was performed inside an inert atmosphere glovebox. Anhydrous FeCl₃ (as-received) and CsCl (dried in air at 400 °C) were loaded into a fused silica crucible. The crucible was placed inside a fused silica ampoule, transferred to a vacuum line without exposure to air, purged with argon and sealed with approximately $\frac{1}{4}$ atm argon. The reagents were heated to 700 °C at a rate of 60 °C/h and held for 8 h prior to quenching in an ice-water bath. A boule was easily extracted from the silica ampoule, which was then ground and the powder was sealed with argon gas in a manner as done for the initial reaction. The powder was annealed at 250 °C for 10 d, with an external thermocouple utilized to verify the annealing temperature; the furnace was turned off and the sample was allowed to cool naturally.

II Single crystal growth

Single crystals of Cs₃Fe₂Cl₉ were synthesized through a solvothermal route with a variety of concentrations and temperature profiles tested. One consistent feature of the growths was that a small excess of CsCl was found to improve the size and apparent quality of the crystals obtained. In the best growth conditions 1 mmol of anhydrous FeCl₃ and 2 mmol anhydrous CsCl were weighed out in a glovebox under argon and removed in a small, closed vial. This vial was then emptied into a 23 mL Parr A280AC PTFE liner in air with concentrated HCl added immediately after, with additional HCl being used as a wash of the vial to ensure all material was transferred, with a total of 5 mL concentrated HCl being added to PTFE liner, which was subsequently sealed in a Parr 4749 general purpose acid digestion vessel. The vessel was placed in an oven which was ramped to 220 °C over one hour and held at this temperature for 6 hours to allow the vessel to come to temperature and allow the reactants to be fully dissolved in the concentrated HCl. The vessel was then cooled to 30 °C over 48 hours before the oven was turned off, with slower cooling rates not providing significant improvements in crystal size.

Due to the pressures involved and the porosity of PTFE the exterior of the liner was coated in a water-soluble green chloride salt due to the corrosion of the acid digestion vessel which was removed before the PTFE liners were opened. The liners were opened in air, with the resulting crystals metastable in the mother liquor, though undergoing a process of dissolution and recrystallization if stored under the mother liquor long-term. To avoid this, and because the crystals were highly air-sensitive they were immediately transferred into concentrated HCl under which they could be safely stored and additionally helped to remove the mother liquor from the crystals. The crystals stored under concentrated HCl were transferred to a glove bag which was filled with argon to create an air free environment safe from corrosion where the crystals were removed from the concentrated HCl, dried via filtration and transferred into a vial under argon which, free of HCl, was then transferred into a helium filled glovebox for long term storage. Crystals obtained through this method were up to 25 mg, though the largest unambiguously single crystal was 7.5 mg which was selected for neutron diffraction experiments on CORELLI and WAND².

III Magnetization measurements

Magnetization measurements were performed on a single crystal weighing 4.5(1) mg using a Quantum Design Magnetic Property Measurement System MPMS-3. Due to the air sensitivity of the crystal it was coated in and secured to a quartz post with silicone vacuum grease. Magnetization measurements revealed similar phase transitions as previously reported Ishii et al. (2021) and were used to plan the neutron diffraction experiment.

IV Neutron diffraction experiments on POWGEN

Powder neutron diffraction experiments were performed on POWGEN Huq et al. (2019) at the Spallation Neutron Source (SNS) of the Oak Ridge National Laboratory (ORNL). About 5 g powder of Cs₃Fe₂Cl₉ was filled into an air-tight vanadium can in a helium filled glovebox. An orange cryostat was utilized to reach a base temperature of 2 K. Data reduction was performed using the MANTID software Arnold et al. (2014).

Figure S1 summarizes the refinement result of the neutron diffraction pattern collected at $T=15$ K. No secondary reflections are observed in the diffraction pattern, which confirms the phase purity of our sample. Refined crystal structure parameters are listed in Table S1. Due to the lacking of neutron diffraction data in the high- $Q$ region, we cannot reliable fit the thermal parameters. Therefore, a uniform thermal factor of $B_{\textrm{iso}}=0.1$ Å² is assumed for all atoms. Figure S2 plots the temperature evolution of the diffraction pattern. Below $T_{\textrm{N}}\sim 5.5$ K, magnetic Bragg peaks belonging to $\bm{q}=(\frac{1}{2},0,0)$ are observed. This transition temperature is consistent with that observed in magnetic susceptibility.

The refinement result of the neutron diffraction pattern collected at $T=1.6$ K is shown in Fig. 1 of the main text. Table S2 summarizes the refined crystal structure parameters. The refined magnitude of the ordered moment is 4.23(6) $\mu_{\textrm{B}}$ .

Table S1: Refined crystal structure parameters for Cs₃Fe₂Cl₉ at 15 K. The space group is

P6_{3}/mmc

with lattice constants

a=7.1784(1)

and

c=17.6374(3)

Å. A uniform thermal parameter of

B_{\textrm{iso}}=0.1

Å² is assumed for all atoms. The goodness-of-fit parameters are

R_{p}=18.6\%

and

R_{wp}=10.2\%

Atom	$x$	$y$	$z$	site
Cs1	0	0	1/4	$2b$
Cs2	1/3	2/3	0.0826(2)	$4f$
Fe	1/3	2/3	0.8458(1)	$4f$
Cl1	0.519(4)	0.0330(3)	1/4	$12j$
Cl2	0.8248(2)	0.6496(4)	0.0892(1)	$12k$

Table S2: Refined crystal structure parameters for Cs₃Fe₂Cl₉ at 1.6 K. The space group is

P6_{3}/mmc

with lattice constants

a=7.1784(1)

and

c=17.6369(2)

Å. A uniform thermal parameter of

B_{\textrm{iso}}=0.1

Å² is assumed for all atoms. The goodness-of-fit parameters are

R_{p}=15.2\%

and

R_{wp}=11.3\%

Atom	$x$	$y$	$z$	site
Cs1	0	0	1/4	$2b$
Cs2	1/3	2/3	0.0829(2)	$4f$
Fe	1/3	2/3	0.8459(1)	$4f$
Cl1	0.510(4)	0.0325(3)	1/4	$12j$
Cl2	0.8247(2)	0.6495(4)	0.0891(1)	$12k$

V Inelastic neutron scattering experiments on CNCS

Inelastic neutron scattering (INS) experiments on powder sample of Cs₃Fe₂Cl₉ were performed on CNCS Ehlers et al. (2011) at the SNS of the ORNL. About 5 g powder was sealed in an aluminum can in a helium filled glovebox. An orange cryostat was utilized to reach a base temperature of 5 K. Measurements were taken with incident neutron energies of $E_{i}=6.59$ , 3.32, and 1.0 meV in the high flux chopper configuration at $T=2$ and 10 K. For each measuring condition, data were collected on an empty can and subtracted as background. Data reduction and projection were performed using the Mslice program in DAVE Azuah et al. (2009).

Figure S3 presents the INS spectra for Cs₃Fe₂Cl₉ measured with an incident neutron nergy of $E_{i}=1.0$ meV at $T=2$ K. An excitation gap of $\sim 0.2$ meV is resolved, which is consistent with the existence of uniaxial SIA in Cs₃Fe₂Cl₉.

VI Diffuse neutron scattering experiments on CORELLI

Single crystal diffuse neutron scattering experiments were performed on CORELLI Ye et al. (2018) at the SNS of the ORNL. A crystal (mass $\sim 7.5$ mg) was aligned with the $c$ axis vertical. A closed cycle refrigerator (CCR) was employed to reach temperatures $T$ down to 5 K. Data were acquired by rotating the sample in 1.5^∘ steps, covering a total range of 360^∘. The counting time at each rotation angle was approximately 1.5 mins with the source operating at a power of 1.4 MW. Data reduction and projection were performed using the MANTID software Arnold et al. (2014). Data collected at $T=20$ K were subtracted as background. Theoretical diffuse neutron scattering patterns were calculated using the self-consistent Gaussian approximation (SCGA) method Conlon and Chalker (2010) as implemented in JuliaSCGA Jul ; Gao et al. (2022).

Figure S4 compares the diffuse scattering patterns calculated using the SCGA method for the fitted $J_{123}$ - $D_{z}$ model at $T=5$ and 6 K. As explained in the main text, the SCGA method is a mean-field-based theory and may underestimate the critical scattering in systems with relatively strong SIA. If we assume $T=6$ K as that in experiments, the calculated scattering pattern will exhibit less contrast at momentum transfers over and out of the spiral surface. Therefore, a reduced temperature of $T=5$ K is assumed in Fig. 3 of the main text to compensate for the underestimated critical scattering. On the other hand, the comparison between the calculations at $T=5$ and 6 K also reveals that the scattering intensity between the integer- $l$ planes is reduced at lower temperatures. This temperature dependence confirms thermal excitations as the origin of the interplanar scattering.

VII Effects of the uniaxial anisotropy on the SSL

Classical Monte Carlo simulations for the fitted $J_{123}$ - $D_{z}$ model were performed to study the effects of the SIA on the SSL in Cs₃Fe₂Cl₉. As a reference, the diffuse scattering pattern at $T=6$ K calculated by the SCGA method is shown in Fig. S5. Through classical Monte Carlo simulations that incorporate critical scattering near $T_{N}$ , Fig. S5(b) summarizes the temperature evolution of the calculated scattering intensity along ( $h+\frac{1}{2}$ , $h+\frac{1}{2}$ , 0), where the integration range is outlined by dashed rectangle in Fig. S5(b). For the fitted $J_{123}$ - $D_{z}$ model, $T_{N}$ is found to be $\sim 5.5$ K. The separated peaks at $T/T_{N}>1$ correspond to thermally stabilized spiral correlations in the SSL state, while their gradually reduced separation with decreasing $T$ reveals the second-order character of the transition into the LRO. Figure S5(d) compares the in-plane ( $x$ ) and out-of-plane ( $z$ ) components of the calculated spin correlations as a function of $T$ . Due to the existence of uniaxial SIA, the evolution of $\langle S_{x}\cdot S_{x}\rangle$ is nonmonotonic, and the ratio between the in-plane and out-of-plane components exhibits an Arrhenius-like behavior in a wide temperature regime.

VIII Neutron diffraction experiments on WAND²

Single crystal neutron diffraction measurements were performed using the WAND² diffractometer Frontzek et al. (2018) at the High Flux Isotope Reactor HFIR with $\lambda=1.486$ Å. A crystal with mass of $\sim 7.5$ mg was aligned with the $c$ axis vertical, and then sealed in an aluminum can with helium exchange gas. A vertical field cryomagnet was used, providing a base temperature of 1.5 K and a maximum field of $H=6$ T. Measurements were performed in phase I at $T$ = 1.5 K and no applied field rotating through 180° in 0.1° steps over 20 h. Measurements in phase II were made at $T$ = 1.5 K and $H=3.2$ T by rotating through 180° in 0.1° steps over 20 h. In phase III at $T=1.5$ K at $H=4.5$ T the sample was rotated through 90° in 0.1° steps over 10 h. In phase IV, at $T=1.5$ K and $H=5.5$ T the sample was rotated through 90° in 0.1° steps over 10 h. In phase V, at $T=4$ K and $H=5.5$ T the sample was rotated through 90° in 0.1° steps over 10 h. In phase VI, at $T=5.05$ K and $H=5.5$ T the sample was rotated over 90° in 0.1° steps with a total measurement time of 26 h.

Figure S6 compares the neutron diffraction patterns in phases III, IV, and V. As explained in the main text, similar diffraction patterns were observed in these three phases. Especially, their magnetic propagation vector stays at the same $\bm{q}_{\textrm{III}}=\frac{5}{6}\times(\frac{1}{3},\frac{1}{3},0)$ , with weak $2\bm{q}_{\textrm{III}}$ higher harmonics observed along the Brillouin zone boundary. Despite their similarities, careful comparisons of the diffraction patterns reveal some differences among these three phases. As emphasized by the dashed rectangle in each panel, the intensity of the $2\bm{q}_{\textrm{III}}$ reflection in phase III is weaker than that in phase IV, which suggests higher field-induced magnetization in phase IV. Comparing phases IV and V, the $\bm{q}_{\textrm{III}}$ reflection becomes weaker in phase V while the $2\bm{q}_{\textrm{III}}$ reflection stays almost unchanged. This evolution is consistent with a transition from a spiral structure in phase IV to a spin density wave (SDW) structure in phase V.

IX Refinements of the nuclear and magnetic structures

Figure S7 compares the experimental intensity of the nuclear Bragg reflections observed on WAND² to the calculations assuming the same crystal structure as that listed in Table S2. The goodness-of-fit parameters are $R_{F2}=26.9\%$ and $R_{F2w}=22.4\%$ . The agreement between the observed and calculated intensity confirms the good quality of our crystal sample.

Figure S8 summarizes the fitting results of the magnetic reflections collected on WAND² in phases I (a), III (b), VI (c), and V (d). The measuring conditions are listed in the main text. In phase VI, only two nonequivalent magnetic reflections were obtained for each magnetic domain, which does not allow a reliable analysis of the magnetic structure. For phase I, the collinear magnetic structure as shown in Fig. 1(a) of the main text was assumed in the fits. The fitted magnitude of the ordered magnetic moment is 2.9(2) $\mu\rm{{}_{B}}$ , with goodness-of-fit parameters $R_{F2}=20.3\%$ and $R_{F2w}=24.4\%$ . For phases III-V, a helical magnetic structure was assumed based on our classical Monte Carlo simulations. The magnitudes of the in-plane magnetic moment, $M_{\parallel}$ , and out-of-plane magnetic moment, $M_{\perp}$ , were treated as fitting parameters. For phase III, the fitted magnitudes are $M_{\parallel}=1.9\pm 1.5$ and $M_{\perp}=1.9\pm 0.4$ $\mu\rm{{}_{B}}$ , with goodness-of-fit parameters $R_{F2}=13.0\%$ and $R_{F2w}=14.7\%$ . For phase IV, the fitted magnitudes are $M_{\parallel}=2.1\pm 1.6$ and $M_{\perp}=1.8\pm 0.5$ $\mu\rm{{}_{B}}$ , with $R_{F2}=9.0\%$ and $R_{F2w}=10.1\%$ . For phase V, the fitted magnitudes are $M_{\parallel}=1.8\pm 4.1$ and $M_{\perp}=1.2\pm 0.8$ $\mu\rm{{}_{B}}$ , with $R_{F2}=3.2\%$ and $R_{F2w}=3.7\%$ . Considering the standard deviations, the fitted magnitude of the in-plane moment $M_{\parallel}$ is consistent with our classical Monte Carlo simulations, where phases III and IV exhibit an elongated helical structure, and phase V exhibits a spin density wave structure with only the $S_{z}$ components.

X Classical Monte Carlo simulations

Classical Monte Carlo simulations for the $J_{123}$ - $D_{z}$ model

\mathcal{H}=\sum_{\langle ij\rangle\in n}J_{n}\bm{S}_{i}\cdot\bm{S}_{j}+D_{z}(% S_{z})^{2}

(S1)

and its perturbed derivatives were performed using the SpinMC code that implements the single spin flip Metropolis algorithm bue . Unless alternately specified, a $12\times 12\times 4$ supercell with 2304 spins was employed in simulations. Observables including the heat capacity and magnetization were averaged over $10^{6}$ of measurement sweeps after $3\times 10^{5}$ thermalization sweeps, where each sweep represents 2304 attempted spin flips at randomly selected sites. Magnetic structure factors were calculated through fast Fourier transform using the FFTW package fft . The parallel tempering algorithm was utilized to facilitate thermal equilibrium, which was performed simultaneously over 128 replicas on 128 cores with a geometric series of temperatures between 1 and 10 K. After every 10 Monte Carlo sweeps, a replica exchange was attempted and the successful exchange rates are above $\sim 25\%$ for all the neighboring replicas.

X.1 The $J_{123}$ model

Figure S9 summarizes the simulation results for the $J_{123}$ model with no SIA. The coupling strengths are $J_{1}=-0.238$ , $J_{2}=0.085$ , $J_{3}=0.053$ meV as fitted from the INS experimental data. In zero field, the $J_{123}$ model exhibits a long-range order transition at $T_{N}\sim 3.5$ K, which is lower than that observed in experiments. The phase diagram shown in Fig. S9(a) mainly consists of two phases. Magnetic structure factors shown in Figs. S9(b-d) with $\bm{q}\sim\frac{4}{3}\times(\frac{1}{3},\frac{1}{3},0)$ in low magnetic fields and $\bm{q}=(\frac{1}{2},\frac{1}{2},0)$ in high magnetic fields, which is very different from the experimental observations.

X.2 The $J_{123}+D_{z}$ model

Figure S10 summarizes the simulation results for the $J_{123}+{D_{z}}$ model that considers isotropic exchange interactions up to the third neighbors and a uniaxial SIA. The coupling strengths are $J_{1}=-0.238$ , $J_{2}=0.085$ , $J_{3}=0.053$ , and $D_{z}=-0.037$ meV as fitted from the INS experimental data. In zero field, a magnetic long-range order transition is observed at $T_{N}\sim 5.8$ K, which is close to the experimental observation. Although the phase diagram shown in Fig. S10(a) mainly consists of two phases as that in Fig. S9(a), magnetic structure factors shown in Figs. S10(b-d) reveal that the propagation vector $\bm{q}$ stays at $(\frac{1}{2},\frac{1}{2},0)$ r.l.u. in both phases. This calculation indicates that the experimentally observed $\bm{q}_{\textrm{I}}=(\frac{1}{2},\frac{1}{2},0)$ r.l.u. in phase I is due to the existence of uniaxial SIA, yet the stabilization of the field-induced phases III, IV, V, and VI requires further perturbations in the spin Hamiltonian.

X.3 An overview of the perturbation effects

As the theoretical phase diagram for the $J_{123}+D_{z}$ model cannot reproduce the experimental results, extensive classical Monte Carlo simulations were performed to explore possible perturbative terms that may stabilize the experimentally observed additional phases in intermediate fields. By confining the perturbations to a single-term Heisenberg exchange interaction up to the tenth neighbors, we calculated the phase diagram and magnetic structure factors for each perturbed model. The exchange paths for $J_{4}$ - $J_{10}$ are indicated in Fig. S11.

Through our calculations, it was found that all of the considered perturbations with a strength lower than $\sim 0.01$ meV are able to introduce additional phases in the intermediate field regime, yet their effects on the magnetic structure factors categorize them into two groups. The first group contains antiferromagnetic (AFM) $J_{6}$ and $J_{7}$ , or ferromagnetic (FM) $J_{8}$ , $J_{9}$ , and $J_{10}$ perturbations. However, in zero field, these types of perturbations do not produce the experimentally observed $\bm{q}_{\textrm{I}}=(\frac{1}{2},\frac{1}{2},0)$ in phase I and thus contradicts the experimental observations. The second group contains AFM $J_{4}$ or FM $J_{5}$ perturbations. This group of perturbations stabilize $\bm{q}_{\textrm{I}}=(\frac{1}{2},\frac{1}{2},0)$ in zero field while introducing additional phases similar to that observed in experiments. For both groups of perturbations, a sign reversal eliminates the field-induced phases, leading to phase diagrams similar to that of the unperturbed model shown in Fig. S10(a).

X.4 The $J_{123}+D_{z}+J_{6}$ model

As an illustration for the first group of perturbations, Fig. S10 summarizes the simulation results for the $J_{123}+{D_{z}}+J_{6}$ model that considers a sixth-neighbor AFM perturbation of $J_{6}=0.01$ meV. Comparison of the phase diagrams shown in Fig. S10(a) and Fig. S12(a) reveals that the perturbation of $J_{6}$ causes additional phases in the intermediate field regime between $\sim 4.5$ and 6.5 T. However, as shown in Fig. S12(b-d), the propagation vector of the long-range ordered phase in zero field deviates from $\bm{q}_{\textrm{I}}=(\frac{1}{2},\frac{1}{2},0)$ r.l.u. and thus fails to reproduce the experimental results. Similar effects are observed for the remaining perturbations in the first group, thus excluding them as the main perturbation to the $J_{123}+{D_{z}}$ Hamiltonian.

X.5 The $J_{123}+D_{z}+J_{5}$ model

As discussed in the previous section, AFM $J_{4}$ and FM $J_{5}$ belong to the second group of perturbations and have similar effects on the phase diagram. However, there is one slight difference for these two perturbations: In zero field, AFM $J_{4}$ and FM $J_{5}$ raises and lowers $T_{N}$ , respectively. This observation favors FM $J_{5}$ since the experimental $T_{N}\sim 5.4$ K is lower than the theoretically predicted $T_{N}\sim 5.8$ K for the unperturbed $J_{123}+D_{z}$ model.

Figures S13-S15 summarize the simulation results for the $J_{123}+{D_{z}}+J_{5}$ model that considers a FM perturbation of $J_{5}=-0.008$ meV. This strength of $J_{5}$ is found to best describe the experimental phase diagram on a $12\times 12\times 4$ supercell. By carefully analyzing the field and temperature dependence of the magnetic susceptibility (see Fig. S14) and specific heat (see Fig. S15), we obtain the phase diagram shown in Fig. S13(a). The magnetic structure factors shown in Figs. S13(b-d) reproduce the $\bm{q}_{\textrm{I}}=(\frac{1}{2},\frac{1}{2},0)$ in zero field and $\bm{q}_{\textrm{III}}=\frac{5}{6}\times(\frac{1}{3},\frac{1}{3},0)$ at $H=5$ T. Therefore, we conclude that FM $J_{5}$ is the main perturbative term to the $J_{123}+{D_{z}}$ Hamiltonian, leading to a $J_{123}+{D_{z}}+J_{5}$ model that is employed in the main text.

X.5.1 Spiral-type orders in phases III and IV

Comparisons of the direction-resolved magnetic structure factors and real-space spin configurations unveil the magnetic structures in each of the field-induced phases. Figure S16 compares the $\langle S_{z}S_{z}\rangle$ and $\langle S_{x}S_{x}\rangle$ components of the magnetic structure factors calculated for phases III, IV, V, and VI through the classical Monte Carlo simulations. In phases III and IV, both the $\langle S_{z}S_{z}\rangle$ and $\langle S_{x}S_{x}\rangle$ components exhibit long-range order, suggesting a spiral-type magnetic order that involves both in-plane and out-of-plane spin components. The magnetic structure factor along the $z$ direction is about 10 times higher than that along the $x$ direction, which indicates that the spiral order in phases III and IV is elongated along the $c$ axis due to the uniaxial SIA. Compared to that in phase III, the relatively stronger $2\bm{q}_{\textrm{III}}$ reflections in phase IV suggests stronger squaring-up effects in a higher magnetic field. This field dependence of the $2\bm{q}_{\textrm{III}}$ reflections is consistent with the experimental observations shown in Fig. S6.

X.5.2 SDW orders in phases V and VI

As shown in Fig. S16(f) and (h), the $\langle S_{x}S_{x}\rangle$ component of the magnetic structure factors in phases V and VI exhibit diffuse patterns, while sharp Bragg peaks are observed in the $\langle S_{z}S_{z}\rangle$ component. This observation indicates that magnetic moments in phases V and VI are ordered only along the $c$ axis, leading to sinusoidally modulated SDW orders in these two phases. We also note that comparison of the $\langle S_{z}S_{z}\rangle$ component in phases III (Fig. S16(a)), IV (Fig. S16(c)), and V (Fig. S16(e)) reveals a similar magnetic propagation vector $\bm{q}_{\textrm{III}}$ and its high harmonics $2\bm{q}_{\textrm{III}}$ , which reproduces the experimental observations shown in Fig. S16.

XI Analytical expressions for the field-induced phases

To further verify the magnetic structures of the field-induced phases, we calculate the magnetic structure factors for the following ansatz on a honeycomb lattice to describe the spiral-type orders in phases III and IV (Eqn. (S2)) and the SDW orders in phases V and VI (Eqn. (S3)):

\bm{M}(\bm{r})=\left\{\begin{array}[]{l}M_{\perp}\cos(\bm{q}\cdot\bm{r})\bm{n}% _{1}+M_{\parallel}\sin(\bm{q}\cdot\bm{r})\bm{n}_{2}+M_{z}\bm{n}_{1}\textrm{,% \qquad\qquad\ \ if}\ \bm{r}\in\{r_{1}\}\\ M_{\perp}\cos(\bm{q}\cdot\bm{r}+\phi)\bm{n}_{1}+M_{\parallel}\sin(\bm{q}\cdot% \bm{r}+\phi)\bm{n}_{2}+M_{z}\bm{n}_{1}\textrm{,\quad if}\ \bm{r}\in\{r_{2}\}% \end{array}\right.

(S2)

\bm{M}(\bm{r})=\left\{\begin{array}[]{l}M_{\perp}\cos(\bm{q}\cdot\bm{r})\bm{n}% _{1}+M_{z}\bm{n}_{1}\textrm{,\qquad\ \ if}\ \bm{r}\in\{r_{1}\}\\ M_{\perp}\cos(\bm{q}\cdot\bm{r}+\phi)\bm{n}_{1}+M_{z}\bm{n}_{1}\textrm{,\quad if% }\ \bm{r}\in\{r_{2}\}\end{array}\right.

(S3)

In these expressions, $M_{\perp}$ and $M_{\parallel}$ are the ordered moments perpendicular to and parallel with the $ab$ plane, respectively. $M_{z}$ is the field induced moment along the $c$ axis. $\bm{n}_{1}$ is a unit vector along the $c$ axis, $\bm{n}_{2}$ is defined in a way that the vectors $\bm{q}$ , $\bm{n}_{2}$ , and $\bm{n}_{1}$ form a cartesian coordinate system. An additional phase factor $\phi$ is introduced for spins on the second sublattice ( $\bm{r}\in\{r_{2}\}$ ).

Figure S17 summarizes the magnetic structure factors for each ansatz. $\bm{M}(\bm{r})$ in real space was first calculated on a $30\times 30$ superlattice, renormalized to equal moment size if $M_{z}\neq 0$ , and then Fourier transferred to reciprocal space. The parameters are $\bm{q}=\bm{q}_{\textrm{III}}=\frac{5}{6}\times(\frac{1}{3},\frac{1}{3})$ , $M_{\perp}=0.9$ , $M_{\parallel}=0.1$ , $\phi=\pi$ , and $M_{z}=0.38$ (0.58) for phase III (IV). The parameters for phase V are $\bm{q}=\bm{q_{\textrm{III}}}=\frac{5}{6}\times(\frac{1}{3},\frac{1}{3})$ , $M_{\perp}=1$ , $M_{z}=0.48$ , and $\phi=\pi$ . The parameters for phase VI are $\bm{q}_{\textrm{VI}}=\frac{11}{9}\times(\frac{1}{3},\frac{1}{3},0)$ , $M_{\perp}=0.5$ , $M_{z}=0$ , and $\phi=0$ . By comparing the magnetic structure factors in Fig. S17 with those in Fig. S16, Fig. S6, and Fig. 4 in the main text, it is confirmed that the double peaks along over the Brillouin zone boundaries are due to high harmonics $2\bm{q}_{\textrm{III}}$ reflections induced by magnetic field, thus corroborating the proposed spiral-type and SDW orders in the field-induced phases.

XII Magnetic order in the magnetization plateau phase

A few possible magnetic orders have been proposed for the 1/2-magnetization plateau phase Ishii et al. (2021). Through classical Monte Carlo simulations, the order in the 1/2-magnetization plateau phase of the $J_{123}+{D_{z}}+J_{5}$ model is presented in Fig. S18. Further neutron diffraction experiments in fields above 6 T will be required to verify this magnetic order.

XIII Impacts of the $J_{5}$ perturbation on the spin dynamics

Figure S19 presents the calculated INS spectra and diffuse neutron scattering pattern in the $(h,k,0)$ plane for the $J_{123}+D_{z}+J_{5}$ model with a perturbative $J_{5}=-0.008$ meV. The INS spectra are calculated by the linear spin wave theory. The diffuse neutron scattering pattern is calculated by the SCGA method at $T=5$ K. The calculated results only exhibit marginal difference compared to those of the $J_{123}+D_{z}$ model as presented in the main text.

References

Ishii et al. (2021) Y. Ishii, Y. Narumi, Y. Matsushita, M. Oda, T. Kida, M. Hagiwara, and H. K. Yoshida, Phys. Rev. B 103, 104433 (2021).
Huq et al. (2019) A. Huq, M. Kirkham, P. F. Peterson, J. P. Hodges, P. S. Whitfield, K. Page, T. Hugle, E. B. Iverson, A. Parizzi, and G. Rennich, J. Appl. Cryst. 52, 1189 (2019).
Arnold et al. (2014) O. Arnold, J. Bilheux, J. Borreguero, A. Buts, S. Campbell, L. Chapon, M. Doucet, N. Draper, R. Ferraz Leal, M. Gigg, V. Lynch, A. Markvardsen, D. Mikkelson, R. Mikkelson, R. Miller, K. Palmen, P. Parker, G. Passos, T. Perring, P. Peterson, S. Ren, M. Reuter, A. Savici, J. Taylor, R. Taylor, R. Tolchenov, W. Zhou, and J. Zikovsky, Nucl. Instr. Methods Phys. Res. Sec. A 764, 156 (2014).
Ehlers et al. (2011) G. Ehlers, A. A. Podlesnyak, J. L. Niedziela, E. B. Iverson, and P. E. Sokol, Rev. Sci. Instrum. 82, 085108 (2011), https://pubs.aip.org/aip/rsi/article-pdf/doi/10.1063/1.3626935/15897689/085108_1_online.pdf .
Azuah et al. (2009) R. T. Azuah, L. R. Kneller, Y. Qiu, P. L. W. Tregenna-Piggott, C. M. Brown, J. R. D. Copley, and R. M. Dimeo, J. Res. Natl. Inst. Stand. Technol. 114, 341 (2009).
Ye et al. (2018) F. Ye, Y. Liu, R. Whitfield, R. Osborn, and S. Rosenkranz, J. Appl. Cryst. 51, 315 (2018).
Conlon and Chalker (2010) P. H. Conlon and J. T. Chalker, Phys. Rev. B 81, 224413 (2010).
(8) “JuliaSCGA,” https://github.com/moon-dust/JuliaSCGA.jl.
Gao et al. (2022) S. Gao, M. A. McGuire, Y. Liu, D. L. Abernathy, C. d. Cruz, M. Frontzek, M. B. Stone, and A. D. Christianson, Phys. Rev. Lett. 128, 227201 (2022).
Frontzek et al. (2018) M. D. Frontzek, R. Whitfield, K. M. Andrews, A. B. Jones, M. Bobrek, K. Vodopivec, B. C. Chakoumakos, and J. A. Fernandez-Baca, Rev. Sci. Instrum. 89, 092801 (2018).
(11) “SpinMC.jl,” https://juliapackages.com/p/spinmc.
(12) “FFTW.jl,” https://github.com/JuliaMath/FFTW.jl.

Codimension-Two Spiral Spin-Liquid in the Effective Honeycomb-Lattice Compound Cs3Fe2Cl9

Abstract

pacs:

Acknowledgements.

References

I Powder sample synthesis

II Single crystal growth

III Magnetization measurements

IV Neutron diffraction experiments on POWGEN

V Inelastic neutron scattering experiments on CNCS

VI Diffuse neutron scattering experiments on CORELLI

VII Effects of the uniaxial anisotropy on the SSL

VIII Neutron diffraction experiments on WAND2

IX Refinements of the nuclear and magnetic structures

X Classical Monte Carlo simulations

X.1 The J123subscript𝐽123J_{123}italic_J start_POSTSUBSCRIPT 123 end_POSTSUBSCRIPT model

X.2 The J123+Dzsubscript𝐽123subscript𝐷𝑧J_{123}+D_{z}italic_J start_POSTSUBSCRIPT 123 end_POSTSUBSCRIPT + italic_D start_POSTSUBSCRIPT italic_z end_POSTSUBSCRIPT model

X.3 An overview of the perturbation effects

X.4 The J123+Dz+J6subscript𝐽123subscript𝐷𝑧subscript𝐽6J_{123}+D_{z}+J_{6}italic_J start_POSTSUBSCRIPT 123 end_POSTSUBSCRIPT + italic_D start_POSTSUBSCRIPT italic_z end_POSTSUBSCRIPT + italic_J start_POSTSUBSCRIPT 6 end_POSTSUBSCRIPT model

X.5 The J123+Dz+J5subscript𝐽123subscript𝐷𝑧subscript𝐽5J_{123}+D_{z}+J_{5}italic_J start_POSTSUBSCRIPT 123 end_POSTSUBSCRIPT + italic_D start_POSTSUBSCRIPT italic_z end_POSTSUBSCRIPT + italic_J start_POSTSUBSCRIPT 5 end_POSTSUBSCRIPT model

X.5.1 Spiral-type orders in phases III and IV

X.5.2 SDW orders in phases V and VI

XI Analytical expressions for the field-induced phases

XII Magnetic order in the magnetization plateau phase

XIII Impacts of the J5subscript𝐽5J_{5}italic_J start_POSTSUBSCRIPT 5 end_POSTSUBSCRIPT perturbation on the spin dynamics

References

Codimension-Two Spiral Spin-Liquid in the Effective Honeycomb-Lattice Compound Cs₃Fe₂Cl₉

VIII Neutron diffraction experiments on WAND²

X.1 The $J_{123}$ model

X.2 The $J_{123}+D_{z}$ model

X.4 The $J_{123}+D_{z}+J_{6}$ model

X.5 The $J_{123}+D_{z}+J_{5}$ model

XIII Impacts of the $J_{5}$ perturbation on the spin dynamics