# Metamagnetism and crystal-field splitting in pseudohexagonal CeRh 3 Si 2

Andrea Amorese, Dmitry Khalyavin, Kurt Kummer, Nicholas B. Brookes, Clemens Ritter, Oksana Zaharko, Camilla Buhl Larsen, Orest Pavlosiuk, Adam P. Pikul, Dariusz Kaczorowski, et al.

► To cite this version:

Andrea Amorese, Dmitry Khalyavin, Kurt Kummer, Nicholas B. Brookes, Clemens Ritter, et al.  
Metamagnetism and crystal-field splitting in pseudohexagonal CeRh 3 Si 2. Physical Review B, 2022,  
105 (12), pp.125119. 10.1103/PhysRevB.105.125119. hal-05074745

HAL Id: hal-05074745

https://cnrs.hal.science/hal-05074745v1

Submitted on 20 May 2025

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
Metamagnetism and crystal-field splitting in pseudohexagonal CeRh₃Si₂

Andrea Amorese,¹,*, Dmitry Khalyavin,² Kurt Kummer,³ Nicholas B. Brookes,³ Clemens Ritter,⁴ Oksana Zaharko,⁵

Camilla Buhl Larsen,⁵ Orest Pavlosiuk,⁶ Adam P. Pikul,⁶ Dariusz Kaczorowski,⁶ Matthias Gutmann,²

Andrew T. Boothroyd,⁷ Andrea Severing,⁸ and Devashibhai T. Adroja,²,⁹,†

¹Max Planck Institute for Chemical Physics of Solids, Nöthnitzer Str. 40, 01187 Dresden, Germany

²ISIS Facility, Rutherford Appleton Laboratory, Chilton, Didcot Oxon OX11 0QX, United Kingdom

³European Synchrotron Radiation Facility, 71 Avenue des Martyrs, CS40220, F-38043 Grenoble Cedex 9, France

⁴Institut Laue-Langevin, 71 avenue des Martyrs, F-38042 Grenoble Cedex 9, France

⁵Laboratory for Neutron Scattering and Imaging, Paul Scherrer Institut, CH-5232 Villigen PSI, Switzerland

⁶Institute of Low Temperature and Structure Research, Polish Academy of Sciences, ul. Okólna 2, 50-422 Wrocław, Poland

⁷Department of Physics, Clarendon Laboratory, Oxford University, Oxford OX1 3PU, United Kingdom

⁸Institute of Physics II, University of Cologne, Zülpicher Str. 77, D-50937 Cologne, Germany

⁹Highly Correlated Matter Research Group, Physics Department, University of Johannesburg, PO Box 524, Auckland Park 2006, South Africa

<figure><img src="image_2.png" /></figure>

(Received 23 August 2021; revised 25 January 2022; accepted 23 February 2022; published 15 March 2022)

CeRh₃Si₂ has been reported to exhibit metamagnetic transitions below 5 K, a giant crystal field splitting, and anisotropic magnetic properties from single crystal magnetization and heat capacity measurements. Here we report results of neutron and x-ray scattering studies of the magnetic structure and crystal-field excitations to further understand the magnetism of this compound. Inelastic neutron scattering and resonant inelastic x-ray scattering reveal a $J_z = 1/2$ ground state for Ce when considering the crystallographic *a* direction as quantization axis, thus explaining the anisotropy of the static susceptibility. Furthermore, we find a total splitting of 78 meV for the *J* = 5/2 multiplet. The neutron diffraction study in zero field reveals that, on cooling from the paramagnetic state, the system first orders at $T_{N1} = 4.7$ K in a longitudinal spin density wave with ordered Ce moments along the *b* axis (i.e., the [0 1 0] crystal direction) and an incommensurate propagation vector **k** = (0, 0.43, 0). Below the lower-temperature transition $T_{N2} = 4.48$ K, the propagation vector locks to the commensurate value **k** = (0, 0.5, 0), with a so-called lock-in transition. Our neutron diffraction study in applied magnetic field *H* || *b* axis shows a change in the commensurate propagation vector and development of a ferromagnetic component at *H* = 3 kOe, followed by a series of transitions before the fully field-induced ferromagnetic phase is reached at *H* = 7 kOe. This explains the nature of the steps previously reported in field-dependent magnetization measurements. A very similar behavior is also observed for the *H* || [0 1 1] crystal direction.

DOI: 10.1103/PhysRevB.105.125119

## I. INTRODUCTION

In cerium-based intermetallic compounds, the interplay between crystal-field effects, hybridization, magnetic ordering, and Kondo screening is responsible for many electronic and magnetic properties that continue to challenge our understanding. In particular, compounds of the CeT₃X₂ family (T = Rh, Ir; X = B, Si) have short Ce-Ce distances and extremely large crystal-field splittings (which gained them the name *giant* crystal-field compounds), which are accompanied by some very unusual magnetic behaviors. The best known compound

in that series is CeRh₃B₂, which crystallizes in a hexagonal (P6/mmm) structure characterized by chains of Ce ions along the hexagonal *c* axis, with extremely short intrachain distances of 3.09 Å [1]. It orders ferromagnetically (with the cerium moments lying in the hexagonal plane) below $T_C = 115$ K, a record-high Curie temperature among cerium intermetallics, but the saturation magnetic moment of 0.4μB/formula unit is strongly reduced compared to what would be expected for a Ce³⁺ ion in a hexagonal crystal field [2–4]. Tentative explanations of such properties have focused on unusually strong hybridization of Ce 4f electrons with the neighboring Rh 4d or Ce 5d states [3], and/or unusually large crystal-field effect that mixes the $J = \frac{5}{2}$ and $J = \frac{7}{2}$ multiplets [5].

The isoelectronic compound CeIr₃B₂, hexagonal at 473 K, adopts a slightly distorted monoclinic C2/m structure at room temperature [6], and is in the valence fluctuating regime, where strong hybridization causes loss of the localized moments [7].

CeRh₃Si₂ and CeIr₃Si₂, on the other hand, crystallize in an orthorhombic (Imma) structure, which is also a slightly distorted version of that of the hexagonal CeRh₃B₂, if their

<footnote>*and.amorese@gmail.com</footnote>

<footnote>†devashibhai.adroja@stfc.ac.uk</footnote>

<footnote>Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI. Open access publication funded by the Max Planck Society.</footnote>
<figure><img src="image_3.png" /><figcaption>FIG. 1. (a) CeRh<sub>3</sub>Si<sub>2</sub> crystal structure, viewed along the orthorhombic *a* axis, which corresponds to the pseudohexagonal *c*<sub>h</sub> axis and is used as the quantization axis *z*′ in the crystal-field analysis. (b) The angular part of the |z̃<sub>2</sub>, ±1/2⟩ ground state orbital, shown on the Ce site.</figcaption></figure>

crystallographic direction *a* is taken as the (pseudo)hexagonal quantization axis (see Fig. 1). Both compounds order antiferromagnetically at low *T* < 5 K, and both show complex behavior that indicates a rearrangement of the moments below *T*<sub>N</sub> at zero applied magnetic field and a cascade of metamagnetic transitions in applied field.

Given the variety of properties arising from similar crystal structures, the CeT<sub>3</sub>X<sub>2</sub> family is a valuable testing ground to study the interplay between unusually large crystal-field effects, hybridization, and magnetic interactions. In particular, it is important to understand (i) the origin of the giant crystal field and its effects on the resulting local Ce states, (ii) the relevant states and mechanisms responsible for the hybridization, and (iii) ultimately the interactions determining the different magnetic structures.

The first step in that direction is the characterization of the electronic states and magnetically ordered phases of Ce. Motivated by this goal, we present microscopic measurements of the 4f crystal-field transitions and magnetic structures in CeRh<sub>3</sub>Si<sub>2</sub>. Its crystal structure and macroscopic magnetic properties were determined by Pikul *et al.* [8]. It crystallizes in the orthorhombic ErRh<sub>3</sub>Si<sub>2</sub>-type structure with lattice constants *a* = 7.133 Å, *b* = 9.734 Å, and *c* = 5.606 Å at 293 K. According to Cenzual *et al.* [9], the structure can be treated as a deformed hexagonal structure with lattice constants *a*<sub>h</sub> and *c*<sub>h</sub> defined by the relations *a* = *c*<sub>h</sub>, *b* = √3*a*<sub>h</sub>, and *c* = *a*<sub>h</sub>. Figure 1 shows that the Ce ions in CeRh<sub>3</sub>Si<sub>2</sub> are surrounded by a distorted hexagonal environment with the hexagonal axis perpendicular to the *bc* plane, and with zigzag chains of Ce atoms along the orthorhombic *a* axis. The interchain Ce–Ce distances (3.58 Å) are fairly short among the intermetallic Ce compounds [10]. The magnetic properties of the system are strongly anisotropic, with moments confined in the orthorhombic *bc* plane with a small in-plane anisotropy that favors *b* as the easy axis [8]. The 4f electrons are localized and the effective magnetic moment at high *T* is close to that calculated for a free Ce<sup>3+</sup> ion. Analysis of the experimental magnetic susceptibility curves in terms of a crystal-field model suggests splittings as large as 60 meV for the ground state *J* = 5/2 multiplet and 240 meV for the higher *J* = 7/2 multiplet.

In zero field, CeRh<sub>3</sub>Si<sub>2</sub> orders antiferromagnetically at *T*<sub>N<sub>1</sub></sub> = 4.7 K with a subsequent first-order phase transition at *T*<sub>N<sub>2</sub></sub> = 4.48 K, which was ascribed to a rearrangement of the antiferromagnetic structure [8]. In the presence of applied magnetic field, these transitions split and evolve independently with temperature [11], and in addition a broad ferromagnetic-like hump appears in the susceptibility plots below *T*<sub>N<sub>2</sub></sub> = 4.48 K. The complex magnetic behavior manifests itself in the isothermal magnetization curves at *T* < *T*<sub>N<sub>2</sub></sub> with field in the *bc* plane as a sequence of distinct steps, which suggests the presence of a cascade of metamagnetic transitions, with a saturation moment of 1.16*μ<sub>B</sub>* per Ce ion [8].

To fully understand the nature of the magnetic ground state and obtain direct information on the magnetic structure and crystal-field energy levels in CeRh<sub>3</sub>Si<sub>2</sub> we performed neutron and x-ray scattering experiments. Neutron diffraction gives insight into the size and direction of the ordered magnetic moment as a function of temperature and magnetic field, and the combination of inelastic neutron scattering (INS) and high-resolution resonant inelastic x-ray scattering (RIXS) provides a complete picture of the crystal-field split multiplet scheme of the Ce<sup>3+</sup> configuration [12–18]. The results obtained from the present study will allow direct comparison with the existing data on CeRh<sub>3</sub>B<sub>2</sub> (Ref. [19]) and CeIr<sub>3</sub>Si<sub>2</sub> (Refs. [20,21]).

## II. EXPERIMENT

### A. Sample preparation

The polycrystalline samples of CeRh<sub>3</sub>Si<sub>2</sub> and LaRh<sub>3</sub>Si<sub>2</sub> used for inelastic neutron scattering and neutron powder diffraction studies were prepared by conventional arc melting of stoichiometric amounts of the constituents Ce (3N purity, Ames Laboratory), Rh ingot (3N purity, Chempur), and Si chips (6N purity, Chempur). A single crystal of CeRh<sub>3</sub>Si<sub>2</sub> was grown by the Czochralski pulling method employing a tetra-arc furnace under protective ultrapure argon atmosphere for RIXS and single crystal neutron diffraction. The starting polycrystalline melt was prepared as mentioned above. The pulling rate was 10 mm/h and the copper heart's rotation speed was 3 rpm. The final ingot of CeRh<sub>3</sub>Si<sub>2</sub> was about 4 mm in diameter and 40 mm in length. The single crystal as well as the polycrystalline samples were wrapped in Ta foil, sealed in an evacuated silica tube, and annealed at 900 °C for 2 weeks. From this single-crystalline rod, a fragment of about 4 mm length, 3 mm width, and 2 mm thickness was cut for neutron diffraction experiments. The quality of the product was verified by means of x-ray powder-diffraction measurements and single crystal Laue diffraction measurements.

### B. Inelastic neutron scattering (INS)

Inelastic neutron scattering measurements on 6.5 g polycrystalline samples of CeRh<sub>3</sub>Si<sub>2</sub> and LaRh<sub>3</sub>Si<sub>2</sub> were performed on the high flux time-of-flight (TOF) spectrometer MERLIN [22] at the ISIS Facility, UK. The samples were placed in an annular geometry thin-wall Al can with outer diameter 40 mm and height 40 mm. A closed-cycle refrigerator (CCR) with helium exchange gas was used to cool the samples down to 5 K. A sloppy chopper was used to select the incident neutron energies of *E*<sub>i</sub> = 100 meV with 350 Hz
frequency and 600 meV with 600 Hz frequency. The elastic
resolution (FWHM) for $E_i = 100$ meV is $\Delta\hbar\omega = 6.5$ meV
(and $\Delta\hbar\omega = 3.3$ meV at $\hbar\omega = 57.0$ meV energy transfer) and
for $E_i = 600$ meV is $\Delta\hbar\omega = 51$ meV (and $\Delta\hbar\omega = 24.0$ meV
at $\hbar\omega = 350$ meV). We measured a standard vanadium sample
under identical conditions to normalize the sample's data in
absolute units of mb/sr/meV/f.u.

C. Resonant inelastic x-ray scattering (RIXS)

The RIXS experiment was performed at the state-of-the-art
soft x-ray beamline ID32 at the European Synchrotron Radiation
Facility (ESRF) in Grenoble, France [23]. The experiment
was performed at the Ce $M_5$ edge detecting the $3d \rightarrow 4f \rightarrow$
$3d$ process, with resonant energy around 882 eV [15]. The
combined beamline resolution during the experiment varied
between 30 meV and 37 meV and was monitored before and
after each acquisition by measuring the (nonresonant) elastic
signal from conductive carbon adhesive tape placed close to
the sample. The carbon tape elastic signal was also used to
determine the zero energy loss position for each spectrum,
with an accuracy of $\pm 2$ meV. The spectra were acquired
in the single photon centroiding mode and treated with the
RIXSToolBox [24]; the energy step size was set to 6.2 meV,
corresponding to splitting each 15 $\mu$m pixel of the Andor
iKon-L CCD detector in 2.7 points, a value tuned to maximize
the performance of the centroiding algorithm [25,26]. The
scattering geometry, the incident photon polarization, and the
incident photon energies were varied for each acquisition to
gain more information by exploiting the cross-section dependence
of the different RIXS excitations. The incident photon
energy was fixed at the central energy of the Ce $M_5$ edge
($E_c = 882.2$ eV) for most acquisitions, and some spectra (see
the Appendix) were measured with 1 eV higher photon energies,
halfway down the high energy tail. A total of 22 RIXS
spectra were acquired with various experimental settings on
three samples. The samples were postcleaved in the cleaving
chamber ($10^{-8}$ mbar) shortly before being transferred to the
main measuring chamber with a base pressure of $\approx 10^{-9}$ mbar.
18 spectra were acquired at the sample temperature of 20 K,
and four more at 300 K, to search for temperature dependencies
of the features and to search for excitations from
thermally populated low-lying excited states. The full set of
acquired spectra is shown in the Supplemental Material [27].

D. Neutron diffraction (ND)

Neutron powder diffraction measurements in zero field
and neutron single crystal diffraction measurements in applied magnetic field on CeRh₃Si₂ were carried out using
the time-of-flight (TOF) neutron diffractometer WISH [28]
at the ISIS Facility, UK. In order to reduce the absorption
from Rh, an annular vanadium can of 8 mm outer diameter
was used. 6.5 g of powder sample was placed in the annular
can with thickness of 1 mm. The sample was cooled to a
base temperature of 1.5 K using a so-called orange cryostat
under helium exchange gas to thermalize the sample. Data
were collected at several temperatures between 1.5 and 6 K.
Rietveld refinements of the crystal and magnetic structures
were performed with the FULLPROF program [29] on the data

measured in detector banks at average $2\theta$ values of 58°, 90°,
122°, and 154°, each covering 32° of the scattering plane.

For single crystal neutron diffraction measurements on the WISH diffractometer in applied magnetic field, the CeRh₃Si₂ crystal (with dimension of 4 mm × 3 mm × 2 mm and a mass of 90 mg) was mounted on an aluminium holder with the *b* axis vertical and the field was applied parallel to the *b* axis. The single crystal data at 1.5 K were collected with five orientations of the *a* axis of the crystal with respect to the incident beam and then projected in the reciprocal plane to identify the magnetic Bragg peaks. We collected data with applied magnetic fields of 0, 3, 4, 6, and 7 kOe. Further, we performed field-dependent neutron diffraction measurements with magnetic fields applied parallel to the [0 1 1] crystal direction (i.e., field in the *bc* plane) using the thermal neutron single crystal instrument ZEBRA at SINQ, PSI, Switzerland. The experiment was performed with neutron wavelength 2.3 Å in the normal beam geometry. We also performed neutron powder diffraction measurements between 2 and 300 K with wavelength 2.41 Å on the D20 diffractometer at Institut Laue-Langevin, Grenoble, to investigate the temperature dependence of the lattice parameters of CeRh₃Si₂ (Appendix 3).

E. Crystal-field model used for the analysis

The spin-orbit interaction splits the 14-fold degenerate energy levels of Ce³⁺ *f*¹ configuration into the sixfold degenerate Hund's rule ground state ²F₃/₂ and the eightfold degenerate ²F₇/₂ multiplet at about 280 meV above the ground state. The crystal field lifts the remaining degeneracies further by splitting the ground state multiplet into three and the excited multiplet into four Kramers doublets. Since the crystal-field potential must reflect the point symmetry of the Ce ion site, it is convenient to express the crystal-field Hamiltonian as a sum of tensor operators $C_n^m$ which transform in the same way as the (renormalized) spherical harmonics, $C_n^m(\theta, \phi) = \sqrt{\frac{4\pi}{2n+1}} Y_n^m(\theta, \phi)$, with parameters $A_n^m$, that are coefficients of $C_n^{-m} + (-1)^m C_n^{+m}$ and have the same normalization as that used by Wybourne in Ref. [30] (see also Refs. [31,32] for more information about the conversion between different crystal-field conventions). For the orthorhombic $C_{2v}$ point symmetry of Ce in CeRh₃Si₂, nine parameters ($A_2^0$, $A_2^2$, $A_4^0$, $A_4^2$, $A_4^4$, $A_6^0$, $A_6^2$, $A_6^4$, and $A_6^6$) have to be determined experimentally, and the resulting crystal-field wave functions are given by linear combinations of the $|J, J_z\rangle$ components that fulfill $\Delta J_z = 2$ and $\Delta J_z = 4$.

The crystal-field analysis can be simplified considerably when treating the structure as pseudohexagonal, as described above (see Fig. 1). Treating the CeRh₃Si₂ as a hexagonal system implies that the *a* axis, with its pseudosixfold rotational symmetry, is used as the quantization axis. We, therefore, express the new coordinate system as (x', y', z') = (-z, y, x) (where the minus can be neglected thanks to the mirror symmetry about the hexagonal plane) in order to avoid confusion with the usual reference system that uses the *c* axis as quantization axis (e.g., Ref. [8]). The hexagonal crystal field only mixes $|J, J_z\rangle$ basis functions that fulfill $\Delta J_z' = 6$, so that the eigenfunctions are the pure $J_z'$ states $|\frac{5}{2}, \pm\frac{1}{2}\rangle$, $|\frac{5}{2}, \pm\frac{3}{2}\rangle$, and $|\frac{5}{2}, \pm\frac{5}{2}\rangle$ for the lowest $2F_{\frac{5}{2}}$ multiplet and $|\frac{7}{2}, \pm\frac{1}{2}\rangle$, $|\frac{7}{2}, \pm\frac{3}{2}\rangle$,
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>FIG. 2. Inelastic neutron scattering data of CeRh$_3$Si$_2$ (black symbols) and LaRh$_3$Si$_2$ (light blue) measured with 100 meV incident energy at 6 K with (a) small and (b) large momentum transfers $Q$. (c) Temperature dependence of the crystal-field excitation at 57.5 meV. (d) Magnetic scattering from CeRh$_3$Si$_2$ at low $Q$ (panel) and high $Q$ (inset) obtained after subtraction of the phonon contribution estimated from the spectrum of LaRh$_3$Si$_2$. The red lines show the crystal-field calculation in the pseudohexagonal approximation (see below) taking into account the dipole and quadrupole cross sections.</figcaption></figure>

$\alpha |\frac{7}{2}, \pm \frac{5}{2}\rangle + \beta |\frac{7}{2}, \mp \frac{7}{2}\rangle$, and $\beta |\frac{7}{2}, \pm \frac{5}{2}\rangle - \alpha |\frac{7}{2}, \mp \frac{7}{2}\rangle$ for the excited $^2F_{\frac{7}{2}}$ multiplet. The active crystal-field parameters are reduced to $A_2^0, A_4^0, A_6^0$, and $A_6^6$. The $A_6^6$ parameter determines the $\alpha$ parameter of the mixed $J_z$ states and is known to have a minor effect on the splitting scheme and on the measurable properties so that it can often be neglected. This further reduces the number of crystal-field parameters to 3. We used inelastic neutron scattering (INS) and resonant inelastic x-ray scattering (RIXS) to determine the $A_n^m$ parameters experimentally.

### III. RESULTS

#### A. Crystal-field scheme

Figure 2 shows the INS data of CeRh$_3$Si$_2$ and of the nonmagnetic reference sample LaRh$_3$Si$_2$, taken just above the magnetic ordering transition. The spectrum of CeRh$_3$Si$_3$ at small momentum transfers, Fig. 2(a), exhibits one sharp excitation at 57.5 meV, which is not visible in the La data nor at large momentum transfer, Fig. 2(b), and must therefore be attributed to a crystal-field excitation. The narrow line width

($\approx$3.87 meV) is resolution limited and in agreement with the localized character of the $f$ electrons in CeRh$_3$Si$_2$. Figure 2(c) shows the temperature dependence of that crystal-field excitation: its intensity decreases, the excitation broadens, and its energy moves towards 55.5 meV as the temperature is increased up to 300 K. These changes are reversed when cooling the sample back down to 6 K. Figure 2(d) shows the data in Fig. 2(a) after subtracting the phonon and background scattering, both estimated from the nonmagnetic reference compound LaRh$_3$Si$_2$. There is tentative evidence for a magnetic peak near 78 meV [see inset to Fig. 2(d)], and also for a magnetic signal between 250 and 450 meV, where transitions between the $J = 5/2$ and $J = 7/2$ multiplets are expected (see Appendix 4). These weak signals in the INS data have insufficient strength to be interpreted conclusively. This is due to limitations of the INS cross section, as discussed below and in Appendix 4.

RIXS, which has recently proved its ability to overcome the cross-section limitations of INS [17], can measure the crystal-field splitting scheme of both the $J = 5/2$ and the $J = 7/2$ multiplets of Ce$^{3+}$ multiplets [16,17], provided the splittings are comparable or larger than the experimental resolution that
<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 3. Experimental RIXS spectra (black) acquired with different experimental configurations: (a)–(c) at 20 K and (d) at 300 K. Spectra (c) and (d) correspond to the same geometry. Inset in (a): isotropic Ce $M_5$ spectrum, showing the incident photon energy used for RIXS. The red dots indicate the transition intensities resulting from the crystal-field analysis and the full red line shows the calculated spectrum considering the present experimental resolution. The elastic intensity at 0 eV, not obtainable with simulations, was adjusted to fit the spectra.</figcaption></figure>

currently amounts to about 30 meV at state-of-the-art beamlines. A selection of three RIXS spectra acquired at a sample temperature of 20 K with various scattering geometries and photon polarizations are shown in Figs. 3(a)–3(c) (black dots). Excitations from the ground state at zero energy transfer into both multiplets are clearly visible: at first glance one peak at about 60 meV, and more than one excitation around 300 meV. Figure 3(d) shows a spectrum acquired in the same geometry as that in (c), but at the higher temperature of 300 K.

The hexagonal crystal-field potential splits the $^2F_{5\overline{2}}$ multiplet into three doublets, one of which will be the ground state, and the $^2F_{7\overline{2}}$ multiplet, centered at about 300 meV, into four doublets; hence we expect to see two excitations from the ground state into the low energy multiplet and four into the higher multiplet. After fitting Voigt profiles to the 18 RIXS spectra measured at 20 K, we found that the transition energies cluster at 59 meV, 78 meV, 265 meV, 308 meV, and 353 meV (see Appendix 1 and Fig. 8 for details). We could not find the sixth transition, which must either be buried under another close-lying excitation or it must have a next-to-zero cross section. Comparing the spectrum at 300 K with the corresponding one at 20 K [Figs. 3(c) and 3(d)] only minor differences are visible in the positive energy loss side, while an anti-Stokes peak appears at about 59 meV negative energy loss. This confirms that the state at 59 meV is the first excited crystal-field state because, if present, the lower lying states

would have been significantly populated at 300 K, thus leading to additional excitations from the excited states. In RIXS we cannot study the temperature dependence of the excitation energies observed in our INS data due to the experimental accuracy of about ±2 meV on the determination of the peak positions.

Having determined the crystal-field transition energies, we now describe the RIXS spectra with a full multiplet RIXS calculation using the QUANTY code [33], including the hexagonal crystal-field approximation. Details about the calculations (reduction factors, optimization of the crystal-field parameters, etc.) can be found in Appendix 2. The following crystal-field parameters best describe the majority of the RIXS spectra: $A_2^0 = -135$ meV, $A_4^0 = -172$ meV, $A_6^0 = 30$ meV, and $A_6^6 = 0$ meV, with a spin-orbit constant $\zeta_{SO} = 76.56$ meV, obtained by reducing the Hartree-Fock atomic parameter to 88% of its value. The resulting crystal-field energy levels and wave functions are listed in Table I. The ground state is a virtually pure $|J_{z'} = \frac{1}{2}\rangle$, followed by the almost pure states $|J_{z'} = \frac{5}{2}\rangle$ and $|J_{z'} = \frac{3}{2}\rangle$. The mixing between the two $J = \frac{5}{2}$ and $J = \frac{7}{2}$ multiplets was found to be very small. For example, the squared amplitudes of the $|\frac{7}{2}, J_{z'}\rangle$ components in the first three crystal-field wave functions amount to 0.0018 for the ground state, 0.02 for the first excited state, and 0.00015 for the second excited state, respectively (see Table IV in Appendix 2). We recall that the notation $J_{z'}$ indicates that the pseudohexagonal
<table><caption>TABLE I. Calculated crystal-field eigenenergies and eigenfunctions in the hexagonal crystal-field approximation with quantization axis $z' \parallel c_h \parallel a$. The state at 345 meV is in parentheses because it has zero or negligible intensity in the RIXS calculations and is not visible in the experimental spectra.</caption><thead><tr><th>|J, J<sub>z'</sub>|</th><th>Energy (meV)</th></tr></thead><tbody><tr><td>|5/2, ±1/2⟩</td><td>0</td></tr><tr><td>|5/2, ±3/2⟩</td><td>59</td></tr><tr><td>|5/2, ±5/2⟩</td><td>78</td></tr><tr><td>|7/2, ±1/2⟩</td><td>265</td></tr><tr><td>|7/2, ±3/2⟩</td><td>308</td></tr><tr><td>|7/2, ±7/2⟩</td><td>(345)</td></tr><tr><td>|7/2, ±5/2⟩</td><td>353</td></tr></tbody></table>

axis $z' \parallel a$ is the quantization axis used to express the crystal field, while previous literature such as Ref. [8] expressed the crystal-field scheme along $z \parallel c$ (see Fig. 1 and the Supplemental Material [27] for further details). The vertical lines in Figs. 3(a)–3(c) represent the peak position and intensities obtained from the crystal-field fits, and the red lines show the calculated crystal-field excitation spectra convoluted with the Voigt line shape. The expected excitation at 345 meV corresponds to the transition $|5/2, ±1/2⟩ \rightarrow |7/2, ±7/2⟩$, i.e., it is forbidden by the RIXS selection rules $\Delta J_{z'} = 0, \pm 1, \pm 2$ and therefore not observable. The elastic intensity (i.e., the difference between the full red line and the dashed red line), being affected by other parameters not included in the calculation (such as the surface roughness), was tuned to best fit the data. The effect of the surface roughness can be clearly seen as an increase in the elastic intensity in the spectrum in Fig. 3(a), which was acquired on a cleaved surface different from that of the spectra in the other panels. The same crystal-field model used to fit the spectra at 20 K also describes very well the data at 300 K, including the anti-Stokes intensity, further confirming the validity of our analysis; see Fig. 3(d).

The proposed crystal-field scheme (Table I) also agrees well with the INS data. The red line in Fig. 2(d) shows the calculated INS cross section including both dipole and quadrupole contributions [34,35], using the same crystal-field scheme obtained by the RIXS analysis. The calculation reproduces the shape of the spectrum, with only one clearly visible peak at 58 meV, but the calculated intensity had to be scaled down by a factor of 2.5 to be comparable with the measured spectrum. This discrepancy in the intensity is due to sample attenuation effects, caused mainly by the relatively large neutron absorption by Rh.

The state at 78 meV has a negligible neutron cross section at low Q because it corresponds to a $|J = 5/2, J_{z'\parallel a} = \pm 1/2⟩ \rightarrow |J = 5/2, J_{z'\parallel a} = \pm 5/2⟩$ transition for the assumed hexagonal point symmetry, and such a transition is dipole forbidden [36]. Therefore, the small but nonzero signal at energies around 78 meV in the low-Q INS spectrum [see Fig. 2(d)] must be either due to the imperfect subtraction of the nonmagnetic (phonon) background or because the true point symmetry is orthorhombic, and therefore some $J_{z'}$ mixing of the states could connect this excitation to the ground state via the dipole interaction. At the same time, the

low intensity at 78 meV confirms the validity of the hexagonal approximation for the crystal-field analysis, meaning that the orthorhombic $J_{z'}$ mixing, if present, is small. To further compare with the RIXS crystal-field scheme, we show INS data at larger momentum transfer [see the inset to Fig. 2(d)], at which quadrupole transitions are expected to gain intensity in INS. The 78 meV transition appears in the beyond-dipole INS simulations (red line in the inset), while it is just barely detectable in the experimental data due to the relatively large noise and background (see Appendix 4 for more details).

## B. Magnetic structure

The powder neutron diffraction patterns collected above $T_{N1} = 4.7$ K can be satisfactorily fitted in the orthorhombic structural model proposed by Cenzual, Chabot, and Parthe for ErRh$_3$Si$_2$ [9]. Below $T_{N1}$, a set of additional diffraction peaks appears, indicating the onset of a long-range magnetic ordering in agreement with the magnetization and specific heat data [8]. The magnetic peaks can be indexed using an incommensurate propagation vector $k = (0, k_y, 0)$ with $k_y \sim 0.439(3)$ at $T = 4.5$ K. A quantitative magnetic structure refinement revealed that the Ce moments are aligned along the b axis forming amplitude modulated longitudinal spin density waves with ferromagnetic coupling along the a and c directions [Fig. 4(b)]. The structure is associated with the two-dimensional $m\Delta_2(\eta, \eta^*)$ irreducible representation of the paramagnetic Imma space group, where $\eta$ and $\eta^*$ are components of the order parameter in the representation space. In the present case, the order parameter is the long-range ordering of magnetic moments localized on the Ce sites and described by a Fourier series as specified below. The symmetry of the magnetic structure is described by the $Icmm1'(0, 0, \gamma)ss0s$ superspace group with the basis vectors and origin related to the paramagnetic structure as $(0, 0, -1, 0)$, $(1, 0, 0, 0)$, $(0, -1, 0, 0)$, $(0, 0, 0, 1)$, and $(-1/4, 1/4, 1/4, 1/4)$, respectively. Below the second magnetic transition at $T_{N2} = 4.48$ K, the propagation vector locks to the commensurate value $k = (0, 1/2, 0)$ [Fig. 4(a)], giving rise to an equal-moment magnetic structure with up-up-down-down (UUDD) stacking of the ferromagnetic (ac) planes as shown in Fig. 4(c). The moment size refined at $T = 1.5$ K is 1.23(3)$\mu_B$. This value of the moment along the $b \parallel y'$ axis is in good agreement with the value $1.286\mu_B$ calculated for a pure $|5/2, \pm 1/2⟩$ ground state (the expectation value of the magnetic moments along the b axis in the presence of a magnetic or exchange field that aligns the moments along $b \parallel y' \perp z'$ is $3/2g_J$ for a $|J = 5/2, J_{z'\parallel a} = \pm 1/2⟩$ ground state doublet). This value is further reduced to about $1.26\mu_B$ if the small intermultiplet mixing is taken into account. The magnetic space group corresponding to the magnetic structure found is 52.314 $P\alpha_nna$. This space group appears with the basis vectors and origin related to the paramagnetic Imma structure as $(0,2,0)$, $(0,0,1)$, $(1,0,0)$, and $(0, 1/2, 0)$, respectively, in the list of $k = (0, 1/2, 0)$, generated by the ISODISTORT [37] and MAXSYM [38] software. The low temperature transition from the incommensurate to commensurate magnetic order is a so-called lock-in transition associated with free-energy terms activated at specific values of propagation vectors. In particular, in the present case with the propagation vector
<figure><img src="image_12.png" /><figcaption>FIG. 4. (a) Rietveld refinement of the neutron diffraction intensity collected on WISH at 1.5 K from the detector bank with average $2\theta = 58.33^{\circ}$. The cross symbols and solid line (red) represent the experimental and calculated intensities, respectively, and the line below (blue) is the difference between them. Tick marks (green) indicate the positions of Bragg peaks: nuclear (top) and magnetic (bottom). The inset shows the specific *d*-spacing range where strong magnetic Bragg peaks are observed. (b) Incommensurate magnetic structure above the lock-in transition. (c) Commensurate magnetic structure below the lock-in transition. The arrows denote the ordered Ce$^{3+}$ magnetic moment directions.</figcaption></figure>

along the $\Delta$ line of symmetry of the *Imma* space group, the free-energy term $\eta^4 + \eta^{*4}$ is allowed when the $k_y$ component takes the commensurate value 1/2. The invariance of this term can be verified using the matrix operators for the generating symmetry elements of the paramagnetic *Imma* space group, listed in Table II. This type of magnetic phase transition has also been observed in some other Ce-based intermetallic compounds such as CeIrGe$_3$ [13] and CeRhGe$_3$ [39].

To explore the origin of the magnetization steps reported by Pikul *et al.* [8], we performed single-crystal neutron diffraction measurements in magnetic field applied along the *b* axis. The (0, *k*, *l*) reciprocal plane demonstrating the field evolution of magnetic satellites is shown in Fig. 5 (top and middle). In agreement with the powder diffraction data, the

zero-field magnetic structure is characterized by the $\mathbf{k} = (0, 1/2, 0)$ propagation vector. The corresponding UUDD magnetic structure with the unit cell twice as big as the nuclear one (see the bottom of Fig. 5) can be presented as a Fourier series $M_{jt} = \sum_{\mathbf{k}} S_{jt} e^{-2\pi i(\mathbf{k}\cdot\mathbf{t}+\phi_{jt})}$ with special choice of the magnetic phases $\phi_{jt}$ as shown in Table III. The series specifies a magnetic moment, $M_{jt}$, localized on the *j*th Ce ion in a unit cell with the lattice translation $\mathbf{t}$. This is the standard way to describe magnetic ordering in crystalline solids [41] and it is a convenient way to discuss the field-induced magnetic structures. Application of the magnetic field $H = 3$ kOe results in a change of the magnetic satellites (Fig. 5). The new satellites reveal the presence of two commensurate propagation vectors $\mathbf{k}_1 = (0, 2/5, 0)$ and $\mathbf{k}_2 = 2\mathbf{k}_1 = (0, 4/5, 0)$. The corresponding magnetic ordering can be interpreted as a UUDD structure whose Fourier decomposition involves these two propagation vectors along with the $\mathbf{k} = 0$ ferromagnetic component (see Table III). Increasing the field up to $H = 4$ kOe reestablishes the $\mathbf{k} = (0, 1/2, 0)$ propagation vector but in addition the second harmonic $\mathbf{k}_2 = 2\mathbf{k}_1 = (0, 1, 0)$, which violates the *I*-centring lattice condition, is also observed. This experimental result reveals another unusual magnetic structure with UUUD stacking of the ferromagnetic layers (Fig. 5). Further increasing the field results in a reappearance of the magnetic satellites with the $\mathbf{k}_1 = (0, 2/5, 0)$ and $\mathbf{k}_2 = 2\mathbf{k}_1 = (0, 4/5, 0)$ propagation vectors and roughly the same intensity at $H = 6$ kOe. The Fourier analysis indicates that these two components, combined with equal amplitudes and an additional $\mathbf{k} = 0$ ferromagnetic moment, constitute a UUUUD magnetic structure (Table III) with the magnetic unit cell five times bigger than the nuclear one (Fig. 5). Finally, increasing the field up to $H = 7$ kOe results in the disappearance of the magnetic satellites, and only a ferromagnetic contribution to the fundamental nuclear reflections satisfying the *I*-centring condition is observed. The observed sequence of the field-induced transitions indicates progressive flipping of the ferromagnetic layers initially aligned antiparallel to the field direction. The number of layers with the moments polarized up becomes bigger, increasing the magnetization along the field direction and resulting in the magnetization steps [8]. The presence of the four distinct field-induced phases with large magnetic unit cells points to the fact that exchange interactions stabilizing these exotic phases have a long-range nature, far beyond the nearest- and next-nearest-neighbor interactions.

The same sequence of magnetic phases were also observed for the field applied along the [011] direction (see Fig. 6), although the critical fields were slightly higher. The orthorhombic structure is a distorted version of the hexagonal lattice shown in Fig. 1(a) with the *bc* plane of the orthorhombic cell corresponding to the *ab* plane of the hexagonal one.

<table><caption>TABLE II. Matrices of the $m\Delta_2(\eta, \eta^*)$ irreducible representation of the *Imma* space group, associated with the $\mathbf{k} = (0, 1/2, 0)$ propagation vector [40]. *T* is the time-reversal operator.</caption><thead><tr><th>{2z|0, 1/2, 0}</th><th>{2y|0, 1/2, 0}</th><th>{-1|0, 0, 0}</th><th>{1|1/2, 1/2, 1/2}</th><th>{1|0, 1, 0}</th><th>T</th></tr></thead><tbody><tr><td>(0 e<sup>3</sup><sub>2</sub>πi e<sup>-1</sup><sub>2</sub>πi)</td><td>(e<sup>1</sup><sub>2</sub>πi 0 e<sup>-1</sup><sub>2</sub>πi)</td><td>(0 1 0)</td><td>(e<sup>1</sup><sub>2</sub>πi 0 e<sup>-1</sup><sub>2</sub>πi)</td><td>(-1 0 -1)</td><td>(-1 0 -1)</td></tr></tbody></table>
<figure><img src="image_13.png" /><figcaption>FIG. 5. (Top) Portion of the (0, *K*, *L*) reciprocal plane measured at 1.5 K in various magnetic fields applied parallel to the *b* axis. (Middle) 1D cut of the reciprocal plane along the [0, *K*, −1] direction. (Bottom) The corresponding magnetic structures.</figcaption></figure>

The two studied field directions correspond to the [120] and [100] directions of the hexagonal lattice and it is not surprising that the response of the crystal is qualitatively similar in both cases.

## IV. DISCUSSION

Our crystal-field analysis of the RIXS and INS spectra found rather large splitting of the low energy $^2F_{5/2}$ multiplet and a pure $|J = \frac{5}{2}, J_z \parallel a \pm \frac{1}{2}\rangle$ ground state with the pseudohexagonal *a* axis as the quantization direction (see the Supplemental Material [27] for the comparison of the two reference systems: $z' \parallel a$ in the present work and $z \parallel c$ in Ref. [8]). Our INS spectra show a resolution-limited crystal-field excitation (Fig. 2) indicating the absence of hybridization of the 4*f* states, in agreement with previous core-level photoemission data [8]. The small temperature dependence of the crystal-field excitation at 58 meV that has been observed with INS is well explained with the thermal expansion of the lattice, as further discussed in Appendix 3.

The deviation from the hexagonal approximation for the crystal field, which would allow some $J_z'$ mixing in all crystal-field levels, must be small. This can be inferred from (i) the quality of the fits shown in Fig. 3 and in the Supplemental Material [27], (ii) the fact that the 78 meV excitation is barely visible in the INS spectra (Fig. 2), and (iii) the nonobservability of the sixth transition in RIXS spectra (see Appendix 1). Furthermore, the proposed crystal-field scheme is in very good agreement with the magnetic susceptibility $\chi(T)$ data from Pikul et al. [8] [Fig. 7(a)]. It is interesting to note that the orthorhombic crystal-field scheme proposed in that study, once rotated onto the $\{|J, J_z'\rangle\}$ basis, corresponds to a 98.2% pure $|\frac{5}{2}, \pm \frac{1}{2}\rangle$ ground state, meaning that the same almost pure ground state was found independently and without making the hexagonal approximation.

The $|\frac{5}{2}, \pm \frac{1}{2}\rangle$ state is spatially elongated along the quantization axis $z' \parallel a$ as shown in Fig. 1(b) and, in the hexagonal approximation, is symmetric about that axis. Such a ground state explains the strong magnetocrystalline anisotropy that constrains the moments in the *bc* easy magnetization plane and is in good agreement with the value of the in-plane
<table><caption>TABLE III. Fourier decomposition $M_{jt} = \sum_k S_{jk}e^{-2\pi i(k \cdot t + \phi_{jk})}$ of the magnetic structures in zero and applied magnetic fields, where $S_{jk}$ is the amplitude of the Fourier component with propagation vector $k$, localized on the $j$th Ce atom ($j = 1, 2$), $\phi_{jk}$ is the magnetic phase, and $t$ is a lattice translation.</caption><thead><tr><th>Field (kOe)</th><th>Propagation vector</th><th>Ce1(0, 1/4, 0.78644)</th><th>Ce2(0, 3/4, 0.21356)</th></tr></thead><tbody><tr><td>0</td><td>$k_1 = (0, 1/2, 0)$</td><td>$S_{1k_1} = 1.75\mu_B$<br>$\phi_{1k_1} = 1/4$</td><td>$S_{2k_1} = 1.75\mu_B$<br>$\phi_{2k_1} = -1/4$</td></tr><tr><td>3</td><td>$k_1 = (0, 2/5, 0)$<br>$k_2 = 2k_1 = (0, 4/5, 0)$<br>$k_3 = 0$</td><td>$S_{1k_1} = 1.59\mu_B$<br>$\phi_{1k_1} = 0$<br>$S_{1k_2} = -0.61\mu_B$<br>$\phi_{1k_2} = 0$<br>$S_{1k_3} = 0.25\mu_B$<br>$\phi_{1k_3} = 0$</td><td>$S_{2k_1} = 1.59\mu_B$<br>$\phi_{2k_1} = 1/5$<br>$S_{2k_2} = -0.61\mu_B$<br>$\phi_{2k_2} = 2/5$<br>$S_{2k_3} = 0.25\mu_B$<br>$\phi_{2k_3} = 0$</td></tr><tr><td>4</td><td>$k_1 = (0, 1/2, 0)$<br>$k_2 = 2k_1 = (0, 1, 0)$<br>$k_3 = 0$</td><td>$S_{1k_1} = 1.23\mu_B$<br>$\phi_{1k_1} = 0$<br>$S_{1k_2} = -0.62\mu_B$<br>$\phi_{1k_2} = 0$<br>$S_{1k_3} = 0.62\mu_B$<br>$\phi_{1k_3} = 0$</td><td>$S_{2k_1} = 1.23\mu_B$<br>$\phi_{2k_1} = 1/4$<br>$S_{2k_2} = -0.62\mu_B$<br>$\phi_{2k_2} = 1/4$<br>$S_{2k_3} = 0.62\mu_B$<br>$\phi_{2k_3} = 0$</td></tr><tr><td>6</td><td>$k_1 = (0, 2/5, 0)$<br>$k_2 = 2k_1 = (0, 4/5, 0)$<br>$k_3 = 0$</td><td>$S_{1k_1} = -0.98\mu_B$<br>$\phi_{1k_1} = 0$<br>$S_{1k_2} = -0.98\mu_B$<br>$\phi_{1k_2} = 0$<br>$S_{1k_3} = 0.74\mu_B$<br>$\phi_{1k_3} = 0$</td><td>$S_{2k_1} = -0.98\mu_B$<br>$\phi_{2k_1} = 1/5$<br>$S_{2k_2} = -0.98\mu_B$<br>$\phi_{2k_2} = 2/5$<br>$S_{2k_3} = 0.74\mu_B$<br>$\phi_{2k_3} = 0$</td></tr><tr><td>7</td><td>$k_1 = 0$</td><td>$S_{1k_1} = 1.23\mu_B$<br>$\phi_{1k_1} = 0$</td><td>$S_{2k_1} = 1.23\mu_B$<br>$\phi_{2k_1} = 0$</td></tr></tbody></table>

moment of $1.23\mu_B$ measured with neutron diffraction and $1.16\mu_B$ estimated from the bulk properties' studies [8]. The value calculated based on our crystal-field model is $1.26\mu_B$. The good agreement of the measured and crystal-field only moment confirms the absence of Kondo screening and also the absence of frustration that potentially could play a role in a pseudohexagonal lattice.

It is worth noting that the same almost pure $|\frac{5}{2}, \pm\frac{1}{2}\rangle$ ground state is found along the pseudohexagonal quantization axis in

<figure><img src="image_14.png" /><figcaption>FIG. 6. Neutron diffraction intensity maps in the $[0, K, 1]$ reciprocal direction as a function of magnetic field ($H$) applied parallel to the $[011]$ direction of CeRh$_3$Si$_2$. The measurements were made on the Zebra diffractometer at $T = 1.5$ K. It is to be noted that the $[011]$ direction is at $30^\circ$ to the $b$ axis.</figcaption></figure>

CeIr$_3$Si$_2$, which has a similar structure and similar magnetic behavior [20,21], and in the strongly hybridized hexagonal ferromagnet CeRh$_3$B$_2$ [19,42]. The crystal-field scheme resulting from our spectroscopic analysis does not support strong effects due to multiplet intermixing in CeRh$_3$Si$_2$, contrary to the expectation from the crystal-field analysis of the static susceptibility [8].

Our neutron diffraction data revealed the appearance of an incommensurate (amplitude modulated) antiferromagnetic structure below $T_{N_1}$, with moments aligned along the orthorhombic $b$ axis, followed by a so-called lock-in of the propagation vector to $k = (0, 1/2, 0)$ below $T_{N_2}$, with a UUDD stacking of the ferromagnetic ac layers. Upon increasing the magnetic field at $T = 1.5$ K, we observe a series of different magnetic structures, alternating between $k = (0, 1/2, 0)$ and $k_1 = (0, 2/5, 0) + k_2$ with $k_2 = 2k_1 = (0, 4/5, 0)$, with, for each structure, a progressive flip of the layers aligned antiparallel to the field direction. These different structures correspond each to a different step seen in the magnetization measurement, explaining the origin of the stairlike shape up to the highest field measured in our experiment. The ratio between the values of the magnetization calculated on the basis of the disproportion between moments pointing “up” and “down” for the different structures at each field value (0 at 0 kOe, 1/5 at 3 kOe, 2/4 at 4 kOe, 3/5 at 6 kOe, and 1 at 7 kOe) perfectly agrees with the ratio between the steps in the magnetization reported in literature [8] [see Fig. 7(b)]. The crystal structure was found to be rigid across all the transitions.

Various rare-earth based systems have been reported to show similar magnetic behavior, with incommensurate antiferromagnetism in the proximity of the Néel temperature
<figure><img src="image_15.png" /><figcaption>FIG. 7. Comparison between the results of this work and the macroscopic measurements from Ref. [8]. (a) Magnetic susceptibility data and crystal-field calculation based on the results of our analysis. Note: the calculated lines for $H \parallel b$ and $H \parallel c$ are perfectly overlapping in the hexagonal crystal-field approximation. (b) Field-dependent magnetization and the expected steps resulting from the structures shown in Fig. 5. The absolute value of the expected magnetization is rescaled to compensate for the incomplete saturation of the moments at 7 kOe in the macroscopic data.</figcaption></figure>

and/or metamagnetic transitions upon changing field strength [43–53]. These phenomena can originate from strong magnetocrystalline anisotropy, competition between long-range ferromagnetic and antiferromagnetic interactions, energy crossing of crystal-field levels upon increasing fields, locking of the moments in intermediate easy directions, or frustration. Recently, the metamagnetic transitions in the famous metamagnetic compound CeSb have been linked to a reorganization of the electronic structure [54].

Based on our data, we can make some considerations on the mechanisms ruling the magnetism in CeRh$_3$Si$_2$. The crystal-field ground state explains the strong magnetic anisotropy. The large crystal-field splitting, on the other hand, excludes crossing or depopulation of crystal-field states upon cooling or applying a magnetic field to be responsible for the

metamagnetic transitions. The collinearity of the magnetism at all temperatures and fields further excludes the locking of the moments in several intermediate easy directions while rotating towards the field.

Having stacks of collinear ferromagnetic planes, and no in-plane propagation vector, one can reduce the problem to a one dimensional Ising chain (along the normal to the planes) with oscillating exchange interactions between nearest- and next-nearest-neighbor planes, and possibly extending even farther. If limited to the next-nearest neighbors, the model is reduced to the so-called axial next-nearest-neighbor Ising model (ANNNI), which was shown able to generate, depending on the parameters, a series of commensurate and incommensurate phases, resulting in the so-called *devil's staircase* shape of the magnetization and in the incommensurate propagation vector near the paramagnetic transition [55–57].

Interestingly, the ANNNI model yields, for a wide range of its parameters, an up-up-down-down ordered ground state at low temperature [55,57,58], i.e., the same magnetic structure found here for CeRh$_3$Si$_2$. Therefore, the metamagnetism is most likely related to oscillating long range exchange interactions, which in rare-earth based compounds can naturally be explained with the Ruderman-Kittel-Kasuya-Yosida interaction.

## V. CONCLUSION

We have investigated the electronic and magnetic structure of CeRh$_3$Si$_2$. Inelastic neutron scattering and resonant inelastic x-ray scattering allowed us to fully characterize the crystal-field scheme of the Ce$^{3+}$ ions, and measured rather large crystal-field splittings of about 58 meV and 78 meV for the Hund's rule ground state multiplet. The INS and RIXS spectra can be well reproduced by simulations using a hexagonal approximation for crystal-field potential (with quantization axis $z'$ $\parallel a$). The details of the spectra have been explained on the basis of the cross section and its selection rules. The resulting crystal-field scheme with a $J_{z'} = \pm \frac{1}{2}$ ground state explains the strong planar anisotropy previously reported in susceptibility measurement, and provides a good agreement with those data and with the measured value of the saturated moment of about 1.23$\mu_B$ along *b*. The residual in-plane anisotropy must be due to the additional orthorhombicity of the system, not considered in our crystal-field analysis.

Neutron diffraction in zero field revealed an incommensurate antiferromagnetic structure below $T_{N1} = 4.7$ K, followed by a locking of the propagation vector below $T_{N2} = 4.48$ K, resulting in an up-up-down-down alternation of the moments, aligned along the *b* axis. In applied field, we identified four different magnetic structures, each corresponding to a step in the field-dependent magnetization. The size of the moment measured is in agreement with the single-ion crystal-field ground state, further confirming the absence of hybridization and/or frustration.

## ACKNOWLEDGMENTS

We would like to thank O. Stockert for many instructive and useful discussions. D.T.A. would like to thank the
<figure><img src="image_16.png" /><figcaption>FIG. 8. Peak positions obtained from Voigt fits to the 18 RIXS spectra at 20 K; yellow symbolizing the transitions into the ground state multiplet and blue the excited multiplet (see Appendix A). The size of the dots corresponds to the relative peak intensity.</figcaption></figure>

Royal Society of London for provision of a Newton Advanced Fellowship for exchange between the UK and China, and International exchange funding between the UK and Japan. We thank the ISIS Facility for beam time on MERLIN, Proposal No. RB2010784 [59], and WISH, Proposal No. RB2010798 [60], SINQ-PSI for beam time on ZEBRA, Proposal ID No. 202100097, and EPSRC UK for the funding grant Reference No. EP/W00562X/1. A.A. and A.S. acknowledge support from the German Research Foundation (DFG), Project No. 387555779. O.P. was supported by the Foundation for Polish Science (FNP), program START No. 66.2020.

# APPENDIX

## 1. Finding transition energies in RIXS spectra

In order to determine the crystal-field splittings from the RIXS data, the positions of the peaks were determined from a Voigt fit using the *fityk* program for data fitting [61]. Each of the 18 spectra acquired at 20 K (the full set of spectra is shown in the Supplemental Material [27]) was fitted independently using five Voigt peaks for the crystal-field excitation, plus one Gaussian function located at 0 meV to account for the elastic signal. The fit was constrained by assuming the same linewidth for all five peaks, with a Gaussian contribution given by the width of the elastic peak. Two peaks were restricted to energies lower than 150 meV, representing the two excited $^2F_{5/2}$ levels, and the other three peaks were assigned in the 150–450 meV region, representing the three visible excitations in the $^2F_{7/2}$ multiplet.

Figure 8 shows the peak positions obtained from all these fits and their relative heights normalized to the sum of the peak heights in each spectrum. The size of the circles is proportional to the absolute intensity of the peaks to give more importance to the fits of spectra with better statistics. The five inelastic peaks are at about the same positions for all the

<figure><img src="image_17.png" /><figcaption>FIG. 9. Regions of the parameter space in agreement with the crystal-field splitting measured by RIXS. The brown spheres indicate the set of crystal field parameters that give the measured splittings, with an error of ±7 meV, the red spheres are the sets of parameters that respect the further conditions on the relative peak intensities, and the yellow spheres respect the stricter condition of having a maximum error of ±2 meV. The $A_6$ parameter (not shown) has a minor effect. To ease the reading of the values, the projection of the spheres are also shown on the planes in lighter colors. The axis unit is eV.</figcaption></figure>

fitted spectra, their intensity-weighted average positions being 0 meV, 59 meV, 78 meV, 265 meV, 308 meV, and 353 meV. The yellow dots in Fig. 8 refer to transitions into the $^2F_{5/2}$ multiplet and the blue ones to transitions into $^2F_{7/2}$. We can observe that independently on the experimental configuration (geometry, polarization, and incident photon energy) the first inelastic peak at 59 meV is stronger than the peak at 78 meV, while the last peak at 353 meV has always very low intensity. As discussed in Appendix 2, the comparison between these observations and the peak ratios in the calculated RIXS spectra can ease the identification of the crystal-field excitations associated to those peaks.

## 2. Full-multiplet calculation of RIXS process including crystal-field model

We performed full-multiplet calculations using the *Quanty* code by Haverkort [33]. In the calculations, the Slater integrals for the electron-electron interactions in the intermediate $3d^9 4f^2$ RIXS state, as well as the $3d$ and $4f$ spin-orbit parameters, were obtained from the Robert D. Cowan Atomic Structure Code [62]. To best reproduce the peak positions, the $4f$ spin-orbit parameter was reduced to 88% of the value calculated for an isolated ion ($\zeta_{SO} = 0.88 \times 0.087$ eV = 0.07656 eV). The $F_{ff}$ integrals for the intermediate RIXS state were scaled to 55% and the $F_{df}$ and $G_{df}$ integrals to 75% of their calculated value, having been tuned by fitting the isotropic XAS spectrum.

We searched in the phase space ($A_2^0, A_4^0, A_6^0, A_6^6$) for sets of crystal-field parameters which correspond to the observed splittings. The parameters were varied in a fine grid of steps and the regions of the phase space compatible with the crystal-field splittings observed in our data are highlighted in Fig. 9.
<table><caption>TABLE IV. Full crystal-field wave functions and energy levels calculated with the crystal-field parameters given in the text, considering the pseudohexagonal quantization axis $z' \parallel c_h \parallel a$. Note that the $A_6^0$ parameter was not determined and kept zero. In the case of finite $A_6^0$, the $J_{z'} = \pm \frac{5}{2}$ and $J_{z'} = \pm \frac{7}{2}$ components are allowed to mix.</caption><thead><tr><th>$|J, J_{z'}|$</th><th>Energy (meV)</th></tr></thead><tbody><tr><td>0.9991 $|\frac{5}{2}, \pm \frac{1}{2}|$ &plusmn; 0.0426 $|\frac{7}{2}, \pm \frac{1}{2}|$</td><td>0</td></tr><tr><td>0.9899 $|\frac{5}{2}, \pm \frac{3}{2}|$ &plusmn; 0.1415 $|\frac{7}{2}, \pm \frac{3}{2}|$</td><td>59</td></tr><tr><td>0.9999 $|\frac{5}{2}, \pm \frac{5}{2}|$ &mp; 0.0122 $|\frac{7}{2}, \pm \frac{5}{2}|$</td><td>78</td></tr><tr><td>0.9991 $|\frac{7}{2}, \pm \frac{1}{2}|$ &mp; 0.0426 $|\frac{5}{2}, \pm \frac{1}{2}|$</td><td>265</td></tr><tr><td>0.9899 $|\frac{7}{2}, \pm \frac{3}{2}|$ &mp; 0.1415 $|\frac{5}{2}, \pm \frac{3}{2}|$</td><td>308</td></tr><tr><td>$|\frac{7}{2}, \pm \frac{7}{2}|$</td><td>345</td></tr><tr><td>0.9999 $|\frac{7}{2}, \pm \frac{5}{2}|$ &plusmn; 0.0122 $|\frac{5}{2}, \pm \frac{5}{2}|$</td><td>353</td></tr></tbody></table>

In particular, the $A_2^0$, $A_4^0$, $A_6^0$, and $A_6^6$ parameters were each varied in the region [-0.8, +0.8] eV in steps of 0.01, 0.02, 0.03, and 0.05 eV, respectively. The brown spheres in Fig. 9 represent the solutions compatible with a $|J = \frac{5}{2}, J_{z'} \parallel a = \pm \frac{1}{2}|$ ground state as determined by magnetic susceptibility studies [8] (once rotated in the $x', y', z'$ reference system used in this work; see Supplemental Material [27] for more details), and the measured splittings, allowing for a ±7 meV error on the position of the crystal-field levels, in order to be sure to take into account all the sets of parameters that yield similar splittings. This region can be further reduced by noticing that in the RIXS calculations the $|\frac{5}{2}, \pm \frac{1}{2}| \to |\frac{5}{2}, \pm \frac{3}{2}|$ excitation always has higher intensity than the $|\frac{5}{2}, \pm \frac{1}{2}| \to |\frac{5}{2}, \pm \frac{5}{2}|$ one, and that the excitation in the $\alpha |\frac{7}{2}, \pm \frac{5}{2}| + \beta |\frac{7}{2}, \mp \frac{7}{2}|$ and $\beta |\frac{7}{2}, \pm \frac{5}{2}| - \alpha |\frac{7}{2}, \mp \frac{7}{2}|$ states has in general a lower intensity relative to the other excitations in the spectrum. Combining this information with the relative peak intensities of Fig. 8, we can assign the 59 meV peak to the $|\frac{5}{2}, \pm \frac{3}{2}|$ state, the 78 meV peak to the $|\frac{5}{2}, \pm \frac{5}{2}|$, and the 353 meV peak to one of the two mixed $\alpha (\beta) \cdot |\frac{7}{2}, \pm \frac{5}{2}| + (-) \beta (\alpha) \cdot |\frac{7}{2}, \pm \frac{7}{2}|$. The red spheres in Fig. 8 indicate the combination of crystal-field parameters complying with these further restrictions. Finally, the yellow spheres show the sets of crystal-field parameters that fit the spectra, obtained by further restricting the condition on the peak positions to a maximum error of ±2 meV. The $A_6^0$ parameter (not shown in Fig. 9) was found to have a minor effect on the calculated spectra. In particular, the most important effect of the $A_6^0$ parameter is to mix the $J_{z'} = \pm \frac{5}{2}$ and $\mp \frac{7}{2}$ components in the two states with highest energy. The measured crystal-field splittings agree with any value of $A_6^0$ in the [-0.05, +0.05] eV range, and from this we can estimate that the maximum allowed $J_{z'}$ mixing is of the order of 25%.

With the crystal-field scheme found, the mixing between the two multiplets is found to be negligible (see Table IV).

### 3. Temperature dependence of crystal-field excitation

In order to understand the origin of the temperature dependence of the crystal-field peak positions of CeRh$_3$Si$_2$, we have measured the temperature dependence of lattice parameters with neutron diffraction on the D20 diffractometer at the Institut Laue-Langevin (ILL), Grenoble, France. The results for

<figure><img src="image_18.png" /><figcaption>FIG. 10. Lattice parameters and unit cell volume as a function of temperature.</figcaption></figure>

the temperature dependence of the lattice parameters and unit cell volume are presented in Fig. 10. All lattice parameters decrease with decreasing temperature from 300 K. Below 20 K, the a lattice parameter is nearly temperature independent, while the b parameter exhibits a small upturn. A smaller upturn was also observed in the c lattice parameter below 10 K. The unit cell volume decreases with temperature from 300 K and becomes temperature independent below 20 K. Considering there is no obvious anomaly observed in the lattice parameters and unit cell volume with temperature, the temperature dependence of the INS peak position could be associated with gradual changes in the lattice and the volume expansion with increasing temperature. In order to estimate the change in the crystal-field peak position, we have used the point charge model to calculate the crystal-field energy using the lattice parameters at 7 K and 300 K. We used the McPhase program for this calculation. Although the first excited crystal-field level estimated from the point-charge model is very high (412.64 meV at 7 K and 410.72 meV at 300 K), the difference between 7 K and 300 K is 1.92 meV, which is in good agreement with the experimental results.

### 4. Details of inelastic neutron scattering cross section and calculations

In Fig. 11 we report the INS spectrum acquired with an incident neutron energy of $E_i = 600$ meV and a calculation of the spectrum from the crystal-field model presented in this work. The calculation includes all relevant multipoles in the neutron-4f scattering potential [32], and the cross section is averaged over the same range of Q as the data. The calculation predicts four excitations into the $J = 7/2$ multiplet between 250 and 450 meV, which are only barely visible in our INS data due to low counting statistics; see Fig. 11(b). The detectability of this part of the spectrum is strongly suppressed by the magnetic form factor and by the broadening given by the resolution, which at these energy transfers is about 28 meV. However, our calculation shows that, with larger samples and longer counting time, it is in principle possible to observe these $J = 5/2$ to $J = 7/2$ excitations by INS.

As discussed in the main text, the excitation near 78 meV is a quadrupolar transition and, therefore, it is not expected
<figure><img src="image_19.png" /><figcaption>FIG. 11. (a) Simulated crystal-field excitation spectrum including multipole interactions beyond the dipole interaction. The spectrum is averaged over a $Q$ range extending from $Q_{min}$ to $10 \AA^{-1}$, where $Q_{min}$ is determined by the neutron kinematics and the lowest scattering angle in the experiment (about 4°). The integrated intensity of each crystal-field transition is written in red above the peak. (b) Comparison between the experimental data (5 K, filled symbols) from MERLIN and the simulated crystal-field spectrum (solid line).</figcaption></figure>

to be visible in the experimental data at the low $Q$ values ($Q < 4 \AA^{-1}$) of the data in Fig. 2(d) or in the calculation of the dipole cross section using the hexagonal CEF approximation, also shown in Fig. 2(d). For larger momentum transfers, quadrupole excitations gain intensity and the structure factor (or form factor) has a maximum at $Q \approx 8 \AA^{-1}$ for $Ce^{3+}$. To compare with the RIXS results we have plotted the INS data at higher momentum transfers, averaged from $Q = 6$ to $11 \AA^{-1}$, and calculated the neutron scattering cross section going beyond the dipole approximation [Fig. 2(d), inset]. The calculation shows how these excitations with $\Delta J_z' = \pm 2$ are accessible at these values of $Q$ and some extra intensity is

<figure><img src="image_20.png" /><figcaption>FIG. 12. Neutron diffraction scans in the [0, K, 1] reciprocal direction as a function of magnetic field (H) applied parallel to the [011] direction of CeRh$_3$Si$_2$. The measurements were made on the Zebra diffractometer at T = 1.5 K. It is to be noted that the [011] direction is at 30° to the b axis.</figcaption></figure>

indeed detected in the experiment around 78 meV, although it is near the detection threshold due to the noise.

### 5. Gaussian fits of neutron diffraction data

In Fig. 12 are shown the 1D cuts of the 2D neutron diffraction intensity map (Fig. 6 of the main text) used to determine the magnetic propagation vectors.

[1] A. C. Lawson, A. Williams, and J. G. Huber, Structural properties of CeRh$_3$B$_2$, J. Less-Common Met. **136**, 87 (1987).

[2] S. K. Dhar, S. K. Malik, and R. Vijayaraghavan, Strong itinerant magnetism in ternary boride CeRh$_3$B$_2$, J. Phys. C: Solid State Phys. **14**, L321 (1981).

[3] J. W. Allen, M. B. Maple, J.-S. Kang, K. N. Yang, M. S. Torikachvili, Y. Lassailly, W. P. Ellis, B. B. Pate, and I. Lindau, Density-of-states-driven transition from superconductivity to ferromagnetism in Ce(Ru$_{1-x}$Rh$_x$)$_3$B$_2$: Scenario for an exchange-split Kondo resonance, Phys. Rev. B **41**, 9013 (1990).

[4] R. Vijayaraghavan, NMR in some mixed valent alloys, J. Magn. Magn. Mater. **52**, 37 (1985).

[5] A. Galatanu, E. Yamamoto, T. Okubo, M. Yamada, A. Thamizhavel, T. Takeuchi, K. Sugiyama, Y. Inada, and Y. Onuki, On the unusual magnetic behaviour of CeRh$_3$B$_2$, J. Phys.: Condens. Matter **15**, S2187 (2003).

[6] K. Kubota, E. Matsuoka, Y. Funasako, T. Mochida, T. Sakurai, H. Ohta, T. Onimaru, T. Takabatake, and H. Sugawara, Weak ferromagnetism below 41 K and structural transition at 395 K in CeIr$_3$B$_2$ single crystal, J. Phys. Soc. Jpn. **82**, 104715 (2013).

[7] K. N. Yang, M. S. Torikachvili, M. B. Maple, and H. C. Ku, Intermediate valence behavior of ternary cerium and uranium transition metal borides, J. Low Temp. Phys. **56**, 601 (1984).

[8] A. P. Pikul, D. Kaczorowski, Z. Gajek, J. Stepień-Damm, A. Ślebarski, M. Werwiński, and A. Szajek, Giant crystal-electric-

field effect and complex magnetic behavior in single-crystalline CeRh$_3$Si$_2$, Phys. Rev. B **81**, 174408 (2010).

[9] K. Cenzual, B. Chabot, and E. Parthé, ErrRh$_3$Si$_2$ and isotopes with an orthorhombic deformation superstructure of the CeCo$_3$B$_2$ type, Acta Crystallogr. Sect. C **44**, 221 (1988).

[10] J. Miyahara, N. Shirakawa, Y. Setoguchi, M. Tsubota, K. Kuroiwa, and J. Kitagawa, Hill plot focusing on Ce compounds with high magnetic ordering temperatures and consequent study of Ce$_2$AuP$_3$, J. Supercond. Novel Magn. **31**, 3559 (2018).

[11] D. Kaczorowski and T. Komatsubara, Complex magnetic behavior in single-crystalline CeRh$_3$Si$_2$, Phys. B: Condens. Matter **403**, 1362 (2008).

[12] P. Burlet, L. P. Regnault, J. Rossat-Mignod, C. Vettier, and J. Flouquet, Neutron scattering and heavy fermion compounds, J. Magn. Magn. Mater. **129**, 10 (1994).

[13] V. K. Anand, A. D. Hillier, D. T. Adroja, D. D. Khalyavin, P. Manuel, G. Andre, S. Rols, and M. M. Koza, Understanding the magnetism in noncentrosymmetric CeIrGe$_3$: Muon spin relaxation and neutron scattering studies, Phys. Rev. B **97**, 184422 (2018).

[14] A. J. Princep, D. Prabhakaran, A. T. Boothroyd, and D. T. Adroja, Crystal-field states of Pr$^{3+}$ in the candidate quantum spin ice Pr$_2$Sn$_2$O$_7$, Phys. Rev. B **88**, 104421 (2013).

[15] A. Amorese, G. Dellea, M. Fanciulli, S. Seiro, C. Geibel, C. Krellner, I. P. Makarova, L. Braicovich, G. Ghiringhelli,
D. V. Vyalikh, N. B. Brookes, and K. Kummer, 4f excitations in Ce Kondo lattices studied by resonant inelastic x-ray scattering, Phys. Rev. B **93**, 165134 (2016).

[16] A. Amorese, N. Caroca-Canales, S. Seiro, C. Krellner, G. Ghiringhelli, N. B. Brookes, D. V. Vyalikh, C. Geibel, and K. Kummer, Crystal electric field in CeRh<sub>2</sub>Si<sub>2</sub> studied with high-resolution resonant inelastic soft x-ray scattering, Phys. Rev. B **97**, 245130 (2018).

[17] A. Amorese, K. Kummer, N. B. Brookes, O. Stockert, D. T. Adroja, A. M. Strydom, A. Sidorenko, H. Winkler, D. A. Zocco, A. Prokofiev, S. Paschen, M. W. Haverkort, L. H. Tjeng, and A. Severing, Determining the local low-energy excitations in the Kondo semimetal CeRu<sub>4</sub>Sn<sub>6</sub> using resonant inelastic x-ray scattering, Phys. Rev. B **98**, 081116(R) (2018).

[18] A. Amorese, O. Stockert, K. Kummer, N. B. Brookes, D.-J. Kim, Z. Fisk, M. W. Haverkort, P. Thalmeier, L. H. Tjeng, and A. Severing, Resonant inelastic x-ray scattering investigation of the crystal-field splitting of Sm<sup>3+</sup> in SmB<sub>6</sub>, Phys. Rev. B **100**, 241107(R) (2019).

[19] F. Givord, J.-X. Boucherle, A. P. Murani, R. Bewley, R.-M. Galéra, and P. Lejay, Crystal electric field excitations in the cerium compound CeRh<sub>3</sub>B<sub>2</sub> studied by inelastic neutron scattering, J. Phys.: Condens. Matter **19**, 506210 (2007).

[20] K. Shigetoh, A. Ishida, Y. Ayabe, T. Onimaru, K. Umeo, Y. Muro, K. Motoya, M. Sera, and T. Takabatake, Easy-plane magnetocrystalline anisotropy in the multistep metamagnet CeIr<sub>3</sub>Si<sub>2</sub>, Phys. Rev. B **76**, 184429 (2007).

[21] Y. Muro, Y. Ohno, T. Okada, and K. Motoya, Multi-step metamagnetism in CeIr<sub>3</sub>Si<sub>2</sub>, J. Magn. Magn. Mater. **310**, 389 (2007).

[22] R. I. Bewley, R. S. Eccleston, K. A. McEwen, S. M. Hayden, M. T. Dove, S. M. Bennington, J. R. Treadgold, and R. L. S. Coleman, MERLIN, a new high count rate spectrometer at ISIS, Phys. B: Condens. Matter **385–386**, 1029 (2006).

[23] N. B. Brookes, F. Yakhou-Harris, K. Kummer, A. Fondacaro, J. C. Cezar, D. Betto, E. Velez-Fort, A. Amorese, G. Ghiringhelli, L. Braicovich, R. Barrett, G. Berruyer, F. Cianciosi, L. Eybert, P. Marion, P. van der Linden, and L. Zhang, The beamline ID32 at the ESRF for soft x-ray high energy resolution resonant inelastic x-ray scattering and polarization dependent x-ray absorption spectroscopy, Nucl. Instrum. Methods Phys. Res., A **903**, 175 (2018).

[24] K. Kummer, A. Tamborrino, A. Amorese, M. Minola, L. Braicovich, N. B. Brookes, and G. Ghiringhelli, RixsToolBox: software for the analysis of soft x-ray RIXS data acquired with 2D detectors, J. Synchrotron Radiat. **24**, 531 (2017).

[25] A. Amorese, C. Langini, G. Dellea, K. Kummer, N. B. Brookes, L. Braicovich, and G. Ghiringhelli, Enhanced spatial resolution of commercial soft x-ray CCD detectors by single-photon centroid reconstruction, Nucl. Instrum. Methods Phys. Res., A **935**, 222 (2019).

[26] A. Amorese, C. Langini, G. Dellea, K. Kummer, N. B. Brookes, L. Braicovich, and G. Ghiringhelli, Enhanced spatial resolution of commercial soft x-ray CCD detectors by single-photon centroid reconstruction: Technical Note, Nucl. Instrum. Methods Phys. Res., A **935**, 227 (2019).

[27] See Supplemental Material at http://link.aps.org/supplemental/10.1103/PhysRevB.105.125119 for the full set of RIXS spectra acquired and the corresponding crystal field calculations, for more information on the crystal field scheme proposed in this

work, and for comparison with the crystal field scheme and reference system used in previous works.

[28] L. C. Chapon, P. Manuel, P. G. Radaelli, C. Benson, L. Perrott, S. Ansell, N. J. Rhodes, D. Raspino, D. Duxbury, E. Spill, and J. Norris, WISH: The new powder and single crystal magnetic diffractometer on the second target station, Neutron News **22**, 22 (2011).

[29] J. R. Carvajal, Recent advances in magnetic structure determination by neutron powder diffraction, Phys. B: Condens. Matter **192**, 55 (1993).

[30] B. G. Wybourne and W. F. Meggers, Spectroscopic properties of rare earths, Phys. Today **18**(9), 70 (1965).

[31] M. Rotter, McPhase User Manual, http://www.mcphase.de/.

[32] A. T. Boothroyd, *Principles of Neutron Scattering from Condensed Matter* (Oxford University Press, Oxford, 2020).

[33] M. W. Haverkort, M. Zwierzycki, and O. K. Andersen, Multi-plet ligand-field theory using Wannier orbitals, Phys. Rev. B **85**, 165113 (2012).

[34] M. Rotter, Using McPhase to calculate magnetic phase diagrams of rare earth compounds, *Proceedings of the International Conference on Magnetism (ICM 2003)* [J. Magn. Magn. Mater. **272–276**, E481 (2004)].

[35] A. T. Boothroyd, SPECTRE: a program for calculating spectroscopic properties of rare-earth ions in crystals (1990–2014). Available at https://xray.physics.ox.ac.uk/software.htm.

[36] I. Zaliznyak and S. Lee, Magnetic Neutron Scattering, in *Modern Techniques for Characterizing Magnetic Materials*, edited by Y. Zhu (Springer, Heidelberg, 2005).

[37] B. J. Campbell, H. T. Stokes, D. E. Tanner, and D. M. Hatch, *ISODISPLACE*: A web-based tool for exploring structural distortions, J. Appl. Crystallogr. **39**, 607 (2006).

[38] J. Perez-Mato, S. Gallego, E. Tasci, L. Elcoro, G. de la Flor, and M. Aroyo, Symmetry-based computational tools for magnetic crystallography, Annu. Rev. Mater. Res. **45**, 217 (2015).

[39] A. D. Hillier, D. T. Adroja, P. Manuel, V. K. Anand, J. W. Taylor, K. A. McEwen, B. D. Rainford, and M. M. Koza, Muon spin relaxation and neutron scattering investigations of the noncentrosymmetric heavy-fermion antiferromagnet CeRhGe<sub>3</sub>, Phys. Rev. B **85**, 134405 (2012).

[40] M. I. Aroyo, A. Kirov, C. Capillas, J. M. Perez-Mato, and H. Wondratschek, Bilbao Crystallographic Server. II. Representations of crystallographic point groups and space groups, Acta Crystallogr. **A62**, 115 (2006).

[41] Y. A. Izyumov, V. E. Naish, and R. P. Ozerov, *Neutron Diffraction of Magnetic Materials* (Consulting Bureau, New York, 1991).

[42] F. Givord, J.-X. Boucherle, R.-M. Galéra, G. Fillion, and P. Lejay, Ferromagnetism and crystal electric field in the cerium compound CeRh<sub>3</sub>B<sub>2</sub>, J. Phys.: Condens. Matter **19**, 356208 (2007).

[43] A. R. Ball, D. Gignoux, A. P. Murani, and D. Schmitt, Crystal field and complex phase diagram in hexagonal PrGa<sub>2</sub>, Phys. B: Condens. Matter **190**, 214 (1993).

[44] A. R. Ball, D. Gignoux, J. Rodriguez Fernandez, and D. Schmitt, Magnetic properties and complex phase diagram of hexagonal NdGa<sub>2</sub>, J. Magn. Magn. Mater. **137**, 281 (1994).

[45] S. L. Bud'ko, Z. Islam, T. A. Wiener, I. R. Fisher, A. H. Lacerda, and P. C. Canfield, Anisotropy and metamagnetism in the RNi<sub>2</sub>Ge<sub>2</sub> (R = Y, La-Nd, Sm-Lu) series, J. Magn. Magn. Mater. **205**, 53 (1999).
[46] T. Shigeoka, H. Fujii, M. Nishi, Y. Uwatoko, T. Takabatake, I. Oguro, K. Motoya, N. Iwata, and Y. Ito, Metamagnetism in TbNi₂Si₂ single crystal, *J. Phys. Soc. Jpn.* **61**, 4559 (1992).

[47] E. Morosan, S. L. Bud'ko, P. C. Canfield, M. S. Torikachvili, and A. H. Lacerda, Thermodynamic and transport properties of RAgGe (R = Tb-Lu) single crystals, *J. Magn. Magn. Mater.* **277**, 298 (2004).

[48] P. A. Goddard, J. Singleton, A. L. Lima-Sharma, E. Morosan, S. J. Blundell, S. L. Bud'ko, and P. C. Canfield, Separation of energy scales in the kagome antiferromagnet TmAgGe: A magnetic-field-orientation study up to 55 T, *Phys. Rev. B* **75**, 094426 (2007).

[49] M. Abliz, K. Kindo, K. Kadowaki, and H. Takeya, High-field magnetization measurements and crystalline electric-field effect in HoNi₂B₂ C, *J. Phys. Soc. Jpn.* **72**, 2599 (2003).

[50] S. Gabáni, S. Maťaš, P. Priputen, K. Flachbart, K. Siemensmeyer, E. Wulf, A. Evdokimova, and N. Shitsevalova, Magnetic structure and phase diagram of TmB₄, *Acta Phys. Pol. A* **113**, 227 (2008).

[51] J. Rossat-Mignod, P. Burlet, S. Quezel, J. M. Effantin, D. Delacôte, H. Bartholin, O. Vogt, and D. Ravot, Magnetic properties of cerium monopnictides, *J. Magn. Magn. Mater.* **31–34**, 398 (1983).

[52] D. Gignoux and D. Schmitt, Rare earth intermetallics, *J. Magn. Magn. Mater.* **100**, 99 (1991).

[53] M. Date, Field-induced magnetic phase transitions, *J. Magn. Magn. Mater.* **90–91**, 1 (1990).

[54] K. Kuroda, Y. Arai, N. Rezaei, S. Kunisada, S. Sakuragi, M. Alaei, Y. Kinoshita, C. Bareille, R. Noguchi, M. Nakayama, S. Akebi, M. Sakano, K. Kawaguchi, M. Arita, S. Ideta, K. Tanaka, H. Kitazawa, K. Okazaki, M. Tokunaga, Y. Haga et al., Devil's staircase transition of the electronic structures in CeSb, *Nat. Commun.* **11**, 2888 (2020).

[55] P. Bak, The devil's staircase, *Phys. Today* **39**(12), 38 (1986).

[56] P. Bak and R. Bruinsma, One-Dimensional Ising Model and the Complete Devil's Staircase, *Phys. Rev. Lett.* **49**, 249 (1982).

[57] P. Bak, Commensurate phases, incommensurate phases and the devil's staircase, *Rep. Prog. Phys.* **45**, 587 (1982).

[58] T. A. Kaplan, Frustrated classical Heisenberg model in one dimension with nearest-neighbor biquadratic exchange: Exact solution for the ground-state phase diagram, *Phys. Rev. B* **80**, 012407 (2009).

[59] A. Pikul, A. Amorese, A. Severing, D. T. Adroja, and D. Kaczorowski, Inelastic Neutron Scattering Study to Investigate Giant Crystal-Electric-Field and Intermultiplet Splitting in CeRh₃Si₂, STFC ISIS Neutron and Muon Source (2020), doi:10.5286/ISIS.E.RB2010784.

[60] D. T. Adroja, D. Khalyavin, A. Pikul, and D. Kaczorowski, Investigation of Multiple Magnetic Structures of CeRh₃Si₂ Using Neutron Diffraction Study, STFC ISIS Neutron and Muon Source (2020), doi:10.5286/ISIS.E.RB2010798.

[61] M. Wojdyr, *Fityk*: a general-purpose peak fitting program, *J. Appl. Crystallogr.* **43**, 1126 (2010).

[62] R. Cowan, *The Theory of Atomic Structure and Spectra* (University of California, Berkeley, 1981).