Accurate simulations of magnetic excitations in the neutron simulation package McStas

Schack, S.B.; Krighaar, K.M.L.; Lenander, E.Y.; Lefmann, K.

doi:10.1107/S1600576726002414

research papers

JOURNAL OF
APPLIED
CRYSTALLOGRAPHY

ISSN: 1600-5767

Volume 59| Part 3| June 2026| Pages 816-826

https://doi.org/10.1107/S1600576726002414

Open

access

Accurate simulations of magnetic excitations in the neutron simulation package McStas

Silas B. Schack,^a Kristine M. L. Krighaar,^a Emma Y. Lenander ^a,^b and Kim Lefmann ^a ^*

^aNanoscience Center, Niels Bohr Institute, University of Copenhagen, Denmark, and ^bInstitute of Physics, École Polytechnique Fédérale de Lausanne (EPFL), CH-1015 Lausanne, Switzerland
^*Correspondence e-mail: [email protected]

Edited by N. B. Christensen, Technical University of Denmark (Received 3 December 2025; accepted 6 March 2026; online 30 April 2026)

This article is part of a collection of articles related to the International Conference on Neutron Scattering, ICNS2025.

An algorithm for the accurate simulation of neutron scattering from magnetic excitations has been developed and implemented as new component for the neutron ray-tracing software McStas. The component SpinWave_BCO simulates inelastic neutron scattering from ferro-, antiferro- and altermagnetic excitations in a body-centred orthorhombic crystal structure, where the dispersion relation and scattered neutron intensities are derived using linear spin wave theory. Data from a simulated triple-axis spectrometer with an extremely high resolution have been verified by direct comparison with theory and by comparison with data simulated using the package SpinW. The component serves as a proof of concept for the implementation of a more general linear spin wave component in McStas.

Keywords: spin waves; neutron scattering; ray-tracing simulations; McStas; altermagnetism.

1. Introduction

Neutron scattering is an essential tool for investigating the structure and dynamics of condensed matter, being sensitive to both the nuclei and the magnetic moments in a sample (Boothroyd, 2020 ). In a neutron scattering experiment, one generally measures both the neutron energy and momentum transfers:

$[\displaystyle\hbar\omega = E_{\rm i}-E_{\rm f} = {{\hbar^{2}} \over {2m_{\rm n}}}(k_{ \rm i}^{2}-k_{\rm f}^{2}),\eqno(1)]$

$[\displaystyle{\bf q} = {\bf k}_{\rm i}-{\bf k}_{\rm f},\eqno(2)]$

where $[m_{\rm n}]$ is the neutron mass, ℏ is Planck's constant, E is the neutron energy, k is the neutron wavevector, and the indices `i' and `f' denote the initial and final neutron states, respectively.

Since access to neutron scattering beam time is limited, much instrumental design and preparation is carried out by neutron ray-tracing simulations. Here, virtual experiments (Lefmann et al., 2008 ) on digital twins of neutron instruments can provide valuable insight into instrument performance, experiment feasibility, experiment resolution, multiple scattering effects and signal-to-background relations. In addition, much simulation work is performed for the analysis and design of novel neutron instrumentation, e.g. at the upcoming European Spallation Source, ESS (Andersen et al., 2020 ).

One of the most used ray-tracing tools is McStas, a package originating from the 1990s (Lefmann & Nielsen, 1999 ) and continuously upgraded since. In the package, neutron instruments are assembled by positioning tailor-made or system-provided components (Willendrup & Lefmann, 2020 ). Among many existing system components, McStas contains a number of model neutron scattering samples, providing incoherent scattering, small-angle scattering, reflectivity and crystal diffraction (Willendrup et al., 2006 ; Willendrup & Lefmann, 2021 ; Christensen et al., 2026 ). However, the number of inelastic scattering samples in McStas is at present very limited, the only examples being a flat inelastic mode (akin to a crystal electric field level), a single quasielastic Lorentzian line and a one-branch phonon sample. This seriously limits the breadth of inelastic virtual experiments that can be performed in McStas.

A component for magnon scattering, called Magnon_bcc, currently exists in the McStas component library and was written by some of us. However, this component is limited in scope and has never been fully developed or validated. In addition, we note that the simulation package McVine contains simulations of general dispersion relations, although the full algorithm is to our knowledge not published (Lin et al., 2016 ; Lin et al., 2019 ).

In this work, we remedy this situation by presenting the new McStas component SpinWave_BCO, which calculates and simulates the spin waves of a ferromagnet or a two-sublattice antiferromagnet in the body-centred orthorhombic structure. The component has been tested in a simple virtual experiment using a thermal triple-axis spectrometer with an unrealistically good resolution in both $[\hbar\omega]$ and q. The simulated dispersion and intensity match convincingly both with analytical calculations and with the output of the package SpinW (Toth & Lake, 2015 ). We expect this component to be useful as a test of instrument design as well as for detailed analysis of spin wave data. The component can be readily upgraded to accept a more general lattice geometry and interaction pattern.

2. Spin wave theory

We here describe the theory of quantized spin waves (magnons), also known as linear spin wave theory, that we have used as a basis for our simulations. Although this topic has been known for more than half a century, spin wave theory comes in different variations, e.g. how to specify anisotropy constants and infamous factor 2 variations in definitions of the interaction constants. For these reasons, we deem it necessary to document exactly which equations lie behind our simulations. The main focus will be on the antiferromagnetic spin waves in the body-centred orthorhombic structure. However, the component will in addition simulate ferromagnetic spin waves in the same crystal structure, so a short introduction to the ferromagnetic case is also given.

The magnetism in the system of interest is governed by Heisenberg interactions with coupling strength J_ij, uniaxial anisotropy parametrized by the constant D and an applied magnetic field with field strength B. For simplicity, we here align both the easy axis and the field along the crystallographic c axis (the z direction). The field is applied to stabilize the spins (s_i) in the positive z direction. The system is thus described by the Hamiltonian

$[{\cal H} = {{1} \over {2}}\sum_{i,j}J_{ij}{\bf s}_{i}\cdot{\bf s}_{j}+\sum_ {i}D(s^{z}_{i})^{2}-g\mu_{\rm B}B\sum_{i}s^{z}_{i},\eqno(3)]$

where is the electron gyromagnetic ratio and $[\mu_{\rm B} = 5.7784\times 10^{-2}]$ meV T⁻¹ is the Bohr magneton. When $[J_{ij}\ \lt\ 0]$ dominates, the system is said to be ferromagnetic, while $[J_{ij}\ \gt\ 0]$ dominating leads to an antiferromagnetic system.

2.1. Ferromagnetic spin waves

The ferromagnetic (FM) case for $[D\ \lt\ 0]$ is relatively simple, as the eigenstate of the Hamiltonian equals the classical ground state, where all spins are oriented along the +z direction.

The magnetic dynamics of the system are found using linear spin wave theory, where the deviations from the ground state are assumed to be small. The Hamiltonian is written using the spin raising and lowering operators for spins at position $[{\bf m}]$ : $[s_{{\bf m}}^{\pm} = s^{x}_{{\bf m}}\pm is_{{\bf m}}^{y}]$ . This is then reduced to second order by the Holstein–Primakoff transformation (Boothroyd, 2020):

$[\displaystyle s_{{\bf m}}^{+} = \sqrt{2S}\sqrt{1-{{a_{{\bf m}}^{\dagger}a_{ {\bf m}}^{\phantom \dagger}}\over{2S}}}\ a_{{\bf m}},\eqno(4)]$

$[\displaystyle s_{{\bf m}}^{-} = \sqrt{2S}\ a^{\dagger}_{{\bf m}}\sqrt{ 1-{{a^{\dagger}_{{\bf m}}a^{\phantom \dagger}_{{\bf m}}}\over{2S}}},\eqno(5)]$

$[\displaystyle s_{{\bf m}}^{z} = S-a^{\dagger}_{{\bf m}}a_{{\bf m}},\eqno(6)]$

where S (unitless) is the spin value, and $[a_{{\bf m}}]$ and its Hermitian conjugate are bosonic operators satisfying the commutation relation $[[a_{{\bf m}}^{\phantom \dagger},a_{{\bf m^{\prime}}}^{\dagger}] = \delta_{{\bf m},{\bf m^{\prime}}}]$ .

Using a first-order expansion of the Holstein–Primakoff transformations in the Hamiltonian, and further Fourier transforming the bosonic operators, leads directly to the diagonal Hamiltonian

$[{\cal H}^{(\rm FM)} = E_{0}+\sum_{{\bf q^{\prime}}}\hbar\omega_{{\bf q}^{\prime} }^{(\rm FM)}a_{{\bf q^{\prime}}}^{\dagger}a_{{\bf q^{\prime}}}.\eqno(7)]$

Here E₀ is the energy of the ground state and $[a_{{\bf q^{\prime}}}^{\dagger}a_{{\bf q^{\prime}}}]$ counts the number of magnons of wavevector $[{\bf q^{\prime}}]$ and energy $[\hbar\omega_{\bf q^{\prime}}]$ , with $[{\bf q^{\prime}}]$ defined within the first Brillouin zone. The spin wave dispersion relation is given by

$[\displaystyle\hbar\omega_{{\bf q}^{\prime}}^{(\rm FM)} = S[J^{(\rm FM)}({\bf 0})-J^ {(\rm FM)}({\bf q^{\prime}})]-\textstyle \left(S-{{1} \over {2}}\right)D+g\mu_{\rm B}B.\eqno(8)]$

Here we have defined the Fourier transform of the coupling constants as

$[\displaystyle J^{(\rm FM)}({\bf q}^{\prime}) = \sum_{{\bolddelta}}J({\bolddelta})\exp({i{\bf q}^{\prime}\cdot{\bolddelta}}),\eqno(9)]$

where $[{\bolddelta}]$ are vectors connecting one particular spin to its neighbours. We assume couplings j_a, j_b and j_c to the nearest neighbours along each crystallographic axis in our body-centred orthorhombic structure. In addition, we add couplings j to the spins in all body-centred positions. The system is shown in Fig. 1 (left). It is convenient to separate $[J^{(\rm FM)}({\bf q)}]$ into two parts accounting for the interactions with the body-centred spins and spins along the crystallographic axes, respectively. We then obtain

$[\displaystyle J^{(\rm FM)}({\bf q^{\prime}}) = J({\bf q^{\prime}})+J_{1}({\bf q^{\prime}}),\eqno(10)]$

with

$[\displaystyle J({\bf q^{\prime})} = 8j\cos\left({{aq^{\prime}_{x}} \over {2}} \right)\cos\left({{bq^{\prime}_{y}} \over {2}}\right)\cos\left({{cq^{\prime}_{ z}} \over {2}}\right)\eqno(11)]$

and

$[J_{1}({\bf q^{\prime}}) = 2[\, j_{a}\cos(aq^{\prime}_{x})+j_{b}\cos(bq^{\prime}_ {y})+j_{c}\cos(cq^{\prime}_{z})].\eqno(12)]$

Neutron scattering from magnons happens at low temperatures primarily via the creation or annihilation of a single magnon. The partial differential cross section for unpolarized neutron scattering from single magnons in the ferromagnet is (Marshall & Lovesey, 1971 )

$[\eqalignno{& \left( {{{\rm d}^{2}\sigma} \over {{\rm d}\Omega \ {\rm d}E_{\rm f}}}\right)_{\pm} = (\gamma r_{ 0})^{2}\left[\,{{g} \over {2}}F({\bf q})\right]^{2}\exp({-2W}){{k_{\rm f}} \over {k_{ \rm i}}}(1+\hat{q}_{z}^{2}){{(2\pi)^{3}} \over {2v_{0}}}S \cr & \quad \quad\times\sum_{{\bf q^{\prime}},{\boldtau}}\left(n_{ {\bf q^{\prime}},a}+\textstyle{{1} \over {2}}\pm{{1} \over {2}}\right)\delta({\bf q}\mp {\bf q^{\prime}}-{\boldtau})\delta(\hbar\omega\mp\hbar\omega^{(\rm FM)}_{{\bf q^{\prime}}}),&(13)}]$

where + () refers to the creation (annihilation) of a magnon with energy $[\hbar\omega^{(\rm FM)}_{{\bf q^{\prime}}}]$ . In the expression, $[\gamma = 1.9130]$ is the neutron gyromagnetic ratio, $[r_{0} = 2.8179]$ fm is the classical electron radius, $[\hat{q}_{z}]$ is the z component of a unit vector along the neutron scattering vector q, v₀ is the volume of the magnetic unit cell and $[n_{{\bf q^{\prime}}} = 1/\{\exp[{\hbar\omega^{(\rm FM)}_{{\bf q^{\prime}}}/(k_ {\rm b}T)}]-1\}]$ is the thermal occupation number of (bosonic) magnons at temperature T, with k_b the Boltzmann constant. The vectors τ are the reciprocal lattice vectors of a sublattice. Both the magnetic form factor, $[F({\bf q})]$ , and the Debye–Waller factor, $[\exp({-2W})]$ , will be approximated as constants equal to 1.

Figure 1
Magnetic structure and exchange parameters for the ferromagnetic (left) and antiferromagnetic (right) system. The plots are generated using SpinW (Toth & Lake, 2015

2.2. Antiferromagnetic spin waves

In this section we describe the dynamics of spin waves in a classical two-sublattice antiferromagnet (the Néel state). The magnetic unit cell of such an antiferromagnetic (AFM) system is shown in Fig. 1 (right).

As with the ferromagnet, the magnetic dynamics are found using linear spin wave theory. In this case we perform the Holstein–Primakoff transformations for each sublattice, where the $[t_{{\bf n}}]$ operators represent the spin operators of the `down' sublattice rotated by 180° about the x axis:

$[\displaystyle s_{{\bf m}}^{+} = \sqrt{2S}\sqrt{1-{{a_{{\bf m}}^{\dagger}a_{ {\bf m}}^{\phantom \dagger}}\over{2S}}}\ a_{{\bf m}}, \quad t_{{\bf n}}^{+} = \sqrt{2S}\sqrt{1-{{b^{\dagger}_{{\bf n}}b^{\phantom \dagger}_{{\bf n}}}\over{2S}}}\ b_ {{\bf n}}, \eqno(14)]$

$[\displaystyle s_{{\bf m}}^{-} = \sqrt{2S}\ a_{{\bf m}}^{\dagger}\sqrt{ 1-{{a_{{\bf m}}^{\dagger}a_{{\bf m}}^{\phantom \dagger}}\over{2S}}},\quad t_{{\bf n}}^{-} \displaystyle = \sqrt{2S}\ b_{{\bf n}}^{\dagger}\sqrt{1-{{b^{\dagger}_{ {\bf n}}b^{\phantom \dagger}_{{\bf n}}}\over{2S}}},\eqno(15)]$

$[\displaystyle s_{{\bf m}}^{z} = S-a^{\dagger}_{{\bf m}}a_{{\bf m}},\quad t^{z }_{{\bf n}} \displaystyle = S-b_{{\bf n}}^{\dagger}b_{{\bf n}}.\eqno(16)]$

After expanding to first order and Fourier transforming to the operators $[S^{\pm}_{{\bf q^{\prime}}}]$ and $[T^{\pm}_{{\bf q^{\prime}}}]$ , the resulting Hamiltonian can be diagonalized by a Bogoliubov transformation to bosonic operators $[\alpha_{{\bf q^{\prime}}}]$ and $[\beta_{{\bf q^{\prime}}}]$ via

$[\eqalign{ S^{+}_{{\bf q^{\prime}}} & = \displaystyle u_{ {\bf q^{\prime}}}\alpha_{{\bf q^{\prime}}}+v_{{\bf q^{\prime}}}\beta_ {{\bf q^{\prime}}}^{\dagger},\cr \displaystyle T^{-}_{{\bf q^{\prime}}} & = \displaystyle u_{{\bf q^{\prime} }}\beta^{\dagger}_{{\bf q^{\prime}}}+v_{{\bf q^{\prime}}}\alpha_{{\bf q}^{\prime}},}\eqno(17)]$

where $[u_{{\bf q^{\prime}}}]$ and $[v_{{\bf q^{\prime}}}]$ are real coefficients fixed by requiring that all operators satisfy the relevant commutation relations, e.g. $[[S^{+}_{{\bf q}^{\prime}},S^{-}_{{\bf q}^{\prime\prime}}] = 2S\delta_{{\bf q^ {\prime}},{\bf q^{\prime\prime}}}]$ and $[[a_{{\bf m}}^{\phantom\dagger},a_{{\bf m^{\prime}}}^{\dagger}] = \delta_{{\bf m},{\bf m^{\prime}}}]$ . This leads to the final Hamiltonian

$[{\cal H} = E_{0}+\sum_{{\bf q^{\prime}}}\hbar\omega_{{\bf q^{\prime}},0 }\left(\alpha_{{\bf q^{\prime}}}^{\dagger}\alpha_{{\bf q^{\prime}}}+ \textstyle{{1} \over {2}}\right)+\sum_{{\bf q^{\prime}}}\hbar\omega_{{\bf q^{\prime}}, 1}\left(\beta_{{\bf q^{\prime}}}^{\dagger}\beta_{{\bf q^{\prime}}}+\textstyle{{1} \over {2}}\right),\eqno(18)]$

where E₀ is the energy of the approximate ground state. $[\hbar\omega_{{\bf q}^{\prime},a}]$ denote the two spin wave dispersion relations, which read (Marshall & Lovesey, 1971)

$[\displaystyle\hbar\omega^{(\rm AFM)}_{{\bf q^{\prime}},a} = \Omega_{{\bf q ^{\prime}}}\pm g\mu_{\rm B}B.\eqno(19)]$

is the mode index, which refers to the positive and negative signs of the Zeeman term, respectively, and

$[\displaystyle\Omega_{{\bf q^{\prime}}} = \sqrt{\big(S\{J({\bf 0})-[J_{1}({\bf 0})-J_{1}({\bf q^{\prime}})]\}-(2S-1)D\big)^{2}-S^2J({\bf q^{\prime}})^{ 2}},\eqno(20)]$

where the Fourier-transformed coupling constants are again given by equations (11) and (12).

The partial differential cross section for unpolarized neutron scattering from single magnons in an antiferromagnet is (Marshall & Lovesey, 1971)

$[\eqalignno{ &\displaystyle\left({{{\rm d}^{2}\sigma} \over {{\rm d}\Omega \ {\rm d}E_{\rm f}}}\right)_{\pm,a} = (\gamma r _{0})^{2}\left[\, {{g} \over {2}}F({\bf q})\right]^{2}\exp({-2W}){{k_{\rm f}} \over {k_{ \rm i}}}(1+\hat{q}_{z}^{2}){{(2\pi)^{3}} \over {4v_{0}}} \cr & \quad\quad\quad\quad\displaystyle\times\sum_{{\bf q^{\prime}},{\boldtau}}\left(n_{ {\bf q^{\prime}},a}+\textstyle{{1} \over {2}}\pm{{1} \over {2}}\right)\delta({\bf q}\mp {\bf q^{\prime}}-{\boldtau})\delta(\hbar\omega\mp\hbar\omega^{(\rm AFM)}_{{\bf q^{\prime}},a})\cr & \quad\quad\quad\quad\displaystyle\times\big{[}u_{{\bf q^{\prime}}}^{2}+v_{{\bf q^{\prime}}}^ {2}+2u_{{\bf q^{\prime}}}v_{{\bf q^{\prime}}}\cos({\boldrho}\cdot {\boldtau})\big{]},&(21)}]$

The factor in the last line of the equation is denoted the coherence factor and consists of the coefficients of the Bogoliubov transformation. The vector $[{\boldrho}]$ connects nearest-neighbour spins on different sublattices and is given by $[{\boldrho} = \textstyle{{1} \over {2}}(a,b,c)]$ for the body-centred structure. One can show that given the condition for conservation of crystal momentum $[{\bf q} = \pm{\bf q^{\prime}}+{\boldtau}]$ the following identity holds:

$[\displaystyle 2u_{{\bf q^{\prime}}}v_{{\bf q^{\prime}}}\cos({ \boldrho}\cdot{\boldtau}) = 2u_{{\bf q}}v_{{\bf q}},\eqno(22)]$

meaning that the coherence factor can be calculated directly from the neutron scattering vector $[{\bf q}]$ . Using this, the coherence factor is given by

$[u_{{\bf q}}^{2}+v_{{\bf q}}^{2}+2u_{{\bf q}}v_{{\bf q}} = 2S\,{{ S[J({\bf 0})\!-\!J({\bf q})\!-\!J_{1}({\bf 0})\!+\!J_{1}({\bf q})]-(2S\!-\!1)D} \over { \Omega_{{\bf q}}}}.\eqno(23)]$

Since $[J({\bf q^{\prime}})]$ has a period of two Brillouin zones, as can be seen from equation (8), the coherence factor causes the scattered intensity to vary between different Brillouin zones, having the highest value when q′ is close to an antiferromagnetic ordering vector.

3. Implementation

We now describe the calculations implemented in the new component SpinWave_BCO in McStas. This component simulates the spin wave spectrum of the magnetic systems described above, incorporating both the dispersion relations (8) and (19) and the scattered intensity, which will be derived from equations (13) and (21).

McStas simulates neutron scattering experiments using Monte Carlo ray-tracing techniques. Neutrons are treated as classical rays defined by the components of their position, velocity and spin. To economize computing power, each ray is treated as a number of neutrons, with the equivalent neutron number given by a (non-integer) weight factor, p, assigned to each ray (Willendrup & Lefmann, 2020). When a ray interacts with a component, the ray parameters, including its weight, are subject to change. The weight factor is modified by the weight multiplier w, such that the weight factor after the interaction with component j is given by

$[p_{j} = w_{j}p_{j-1}.\eqno(24)]$

Scattering events are simulated using Monte Carlo sampling, where we denote the sampling probability of a given event a by $[f_{\rm MC}^{a}]$ . The weight multiplier w_a is then the ratio between the physical probability of that event P^a, which is found from the cross section, and the sampling probability $[f_{\rm MC}^{a}]$ (Willendrup & Lefmann, 2021):

$[\displaystyle w_{a} = {{P^{a}} \over {f_{\rm MC}^{a}}}.\eqno(25)]$

As an example, to simulate reflection in a neutron guide, one would only simulate reflected neutrons, meaning $[f_{\rm MC}^{\rm ref.} = 1]$ . If e.g. the probability of reflection is $[P^{\rm ref.} = 0.8]$ , we would obtain a weight multiplier of $[w_{\rm ref.} = 0.8]$ .

For the SpinWave_BCO component, the Monte Carlo choices consist of (a) selecting a scattering direction, $[\hat{\bf k}_{\rm f}]$ , and (for the antiferromagnetic system) (b) selecting one spin wave mode out of $[n_{\rm m} = 2]$ possibilities. The former is determined with a uniform distribution within a predefined solid-angle interval, $[\Delta\Omega]$ . When a neutron scatters from an excitation following a dispersion relation $[\hbar\omega_{{\bf q}}]$ , in this case a magnon in the chosen magnetic system, it must obey the kinematic constraint

$[\displaystyle|\hbar\omega| = \hbar\omega_{{\bf q}},\eqno(26)]$

where the spin wave dispersion $[\hbar\omega_{{\bf q^{\prime}}}]$ is now generalized to the neutron wavevector $[{\bf q}]$ by periodicity. Via the kinematic constraint, the neutron fulfils one of the delta functions in equation (13) or (21). Given that the direction of $[{\bf k}_{\rm f}]$ is already chosen, the kinematic constraint can be solved, giving $[n_{\rm f}\geq 1]$ possible values of the final neutron speed $[v_{\rm f}]$ (Squires, 1977 ). These values of $[v_{\rm f}]$ are determined, and a third MC choice (c) selects one of these values. This algorithm is essentially equivalent to MC simulation of phonon scattering, implemented in McStas through the components Phonon_simple and Phonon_PG and described by Davidsen et al. (2026 ). It is here derived that the weight multiplier of the component is given by

$[w = {{n_{\rm f}n_{\rm m}\Delta\Omega l_{\rm max}} \over {Nv_{0}}}{{{\rm d}\sigma} \over {{\rm d}\Omega}},\eqno(27)]$

where $[l_{\rm max}]$ is the path length of that particular neutron ray inside the sample if it does not scatter, Nv₀ is the total volume of the sample and $[{\rm d}\sigma/{\rm d}\Omega]$ is the differential cross section. The differential cross section is defined as

$[{{{\rm d}\sigma} \over {{\rm d}\Omega}} = {{{\rm Intensity\ of\ neutrons\ scattered\ into\ } {\rm d}\Omega} \over {\Psi \ {\rm d}\Omega}},\eqno(28)]$

where Ψ is the flux of the incoming neutron beam. The cross section thus describes the total number of neutrons scattered into a particular direction, regardless of their energy (Boothroyd, 2020).

When searching for solutions to the kinematic constraint, the initial neutron velocity and the direction of the neutrons' final velocity, $[{\bf v}_{\rm i}]$ and $[{\hat{\bf v}}_{\rm f}]$ , respectively, are kept fixed, while the neutron scattering vector and energy transfer are defined by the final neutron speed according to

$[\displaystyle{\bf q}(v_{\rm f}) = {{m_{\rm n}} \over {\hbar}}({\bf v}_{\rm i}-v_{ \rm f}{\hat{\bf v}}_{\rm f}),\eqno(29)]$

$[\displaystyle\hbar\omega(v_{\rm f}) = {{m_{\rm n}} \over {2}}(v_{\rm i}^{2}-v_{\rm f}^{2 }).\eqno(30)]$

This transforms the kinematic constraint into a one-dimensional problem, easily solved by standard numerical techniques; for details of this procedure, see Davidsen et al. (2026). The differential cross section $[{\rm d}\sigma/{\rm d}\Omega]$ is calculated from $[{\rm d}^{2}\sigma/({\rm d}\Omega \ {\rm d} E_{\rm f})]$ [equations (13) and (21) for the ferromagnet and antiferromagnet, respectively] by integrating over $[E_{\rm f}]$ . To do this, we utilize the relation

$[\int\limits^{\infty}_{0}{\rm d}x \ g(x)\delta[f(x)] = \sum_{j = 1}^{n_{\rm f}}{{g(x_{j})} \over { \displaystyle \left|{({{\rm d}f} / {{\rm d}x}})(x_{j})\right|}} = \sum_{j = 1}^{n_{\rm f}}{{g(x_{j})} \over {| {\cal J}(x_{j})|}}.\eqno(31)]$

Here, x_j is the jth (out of $[n_{\rm f}]$ ) positive zero of f(x), and we have defined $[{\cal J}(x_{j}) = |({\rm d}f/{\rm d}x)(x_{j})|]$ . We can now perform the integration of the inelastic cross section, leading to the expression for the differential cross section for the creation (annihilation) of a single magnon with final speed $[v_{{\rm f},j}]$ out of the $[n_{\rm f}]$ possible values. For the ferromagnetic case, this leads to

$[\eqalignno{ \displaystyle\left({{{\rm d}\sigma} \over {{\rm d}\Omega}}\right)_{\pm,j} = & \ (\gamma r_{0})^{2} \left[\, {{g} \over {2}}F({\bf q})\right]^{2}\exp({-2W}){{{v_{{\rm f},j}}} \over {v_{\rm i }}}[1+(\hat{q}_{j}^{z})^{2}]{{N} \over {2}}S \cr & \displaystyle\times\left(n_{{\bf q}_{j}}+\textstyle{{1} \over {2}}\pm{{1} \over {2}}\right) {{m_{\rm n}v_{{\rm f},j}} \over {|{\cal J}(v_{{\rm f},j})|}},&(32)}]$

while for the antiferromagnetic case, the expression for mode a is

$[\eqalignno{ \displaystyle\left({{{\rm d}\sigma} \over {{\rm d}\Omega}}\right)_{\pm, a,j} = &\ (\gamma r_{0})^{2} \left[\,{{g} \over {2}}F({\bf q})\right]^{2}\exp(-2W){{{v_{{\rm f},j}}} \over {v_{\rm i }}}[1+(\hat{q}_{j}^{z})^{2}]{{N} \over {4}} \cr & \displaystyle\times\left(n_{{\bf q}_{j},a}+\textstyle{{1} \over {2}}\pm{{1} \over {2}}\right) (u_{{\bf q}_{j}}^{2}+v_{{\bf q}_{j}}^{2}+2u_{{\bf q}_{j}}v_{{\bf q }_{j}}){{m_{\rm n}v_{{\rm f},j}} \over {|{\cal J}_{a}(v_{{\rm f},j})|}}.\cr & &(33)}]$

In the simulation, $[{\cal J}(v_{\rm f}) = ({{\partial} /{\partial v_{\rm f}}})[|\hbar\omega(v_{\rm f })|-\hbar\omega_{{\bf q}(v_{\rm f})}]]$ is found numerically for the selected value of $[v_{\rm f}]$ .

The input parameters for the component are listed in Table 1.

Table 1
Table of input parameters for the spin wave component SpinWave_BCO

Parameter name	Description
`radius`, `yheight` [m]	Radius and height of cylindrical sample.
`sigma_inc`, `sigma_abs` [barns]	Cross sections for incoherent scattering and absorption.
`target_index`	Relative index of component to focus at with next component being +1.
`target_x`, `target_y`, `target_z` [m]	Coordinates to focus at (alternative to `target_index`).
`focus_xw`, `focus_yh` [m]	Width and height of focusing.
`focus_aw`, `focus_ah` [°]	Direct specification of ΔΩ (alternative to `focus_xw` and `focus_yh`)
`FM`	Flag for type of order. Default is 0 (antiferromagnet), 1 selects ferromagnetic.
`a`, `b`, `c` [Å]	Lattice constants of orthorhombic lattice unit cell.
`S`	Spin of magnetic ions.
`j` [meV]	Magnetic interaction constant to the eight equivalent neighbours in the body-centred positions.
`j_a`, `j_b`, `j_c` [meV]	Magnetic interaction constants to the neighbours along the a, b and c axes, respectively.
`j_110`, `j_110_prime` [meV]	Magnetic interaction constants for altermagnetic splitting when `FM` = 0 (see Section 5). Both are 0 by default.
`D` [meV]	Uniaxial single-ion anisotropy constant.
`B` [T]	Magnetic field strength.
`T` [K]	Temperature.
`mode_input`	If `FM` = 0: Parameter to specify which mode is simulated. Value of 2 (default) simulates both modes. Value of 0 or 1 simulates only the corresponding mode.

The component takes the form of a cylindrical sample. When it is placed without rotation, the a, b and c axes lie along the x, y and z axes of the McStas coordinate system, respectively.

The main functionality of the component is listed in Algorithm 1:

In these calculations, by far the largest computational power is expended in finding the relevant values of $[v_{\rm f}]$ that simultaneously fulfil the dispersion relation and the kinematic constraint.

The component also conducts a number of checks of the calculated expressions to ensure that the input parameters lead to physically meaningful results. An error message is produced if the dispersion becomes negative or imaginary and if the coherence factor becomes negative when simulating an antiferromagnetic system.

These would probably be caused by the input parameters that do not correspond to the chosen ground state (FM or AFM). When simulating a ferromagnetic system, the value of B is forced to be positive (or zero), such that the applied field always stabilizes the ground state. For the antiferromagnetic system, the magnetic field strength is set to zero if it induces a spin-flop transition, i.e. if the dispersion for the low-energy mode becomes negative.

4. Validation

The component was validated by comparing theoretical calculations with simulated data. This includes a validation of both the simulated dispersion relation, by extracting the Hamiltonian parameters from simulated data, and the simulated absolute intensities.

The main focus was on the validation of the antiferromagnet. This was done by modelling the spin wave spectrum of MnF₂. This material has previously been understood as a classical antiferromagnet, with the magnetic Mn²⁺ ion situated in a body-centred tetragonal lattice. The crystal and Hamiltonian parameters for MnF₂ used in the simulations are listed in Table 2(b).

Table 2
Input parameters to the spin wave component used in the validation process

The antiferromagnetic parameters (b) match those modelled for the compound MnF₂; adopted from Okazaki et al. (1964) and Nikotin et al. (1969).

(a)		(b)
Ferromagnet		Antiferromagnet
Parameter	Value	Parameter	Value
	4.873 Å		4.873 Å
c	3.130 Å	c	3.130 Å
S	2.5	S	2.5
j	−0.304 meV	j	0.304 meV
j_c	−0.056 meV	j_c	−0.056 meV
$[j_{a} = j_{b}]$	−0.008 meV	$[j_{a} = j_{b}]$	0.008 meV
D	−0.023 meV	D	−0.023 meV

We employ a virtual thermal triple-axis spectrometer (TAS), operating with a fixed final energy of $[E_{\rm f} = 14.7]$ meV. The instrument was previously used to test the Phonon_PG component and is described by Davidsen et al. (2026). The simulated instrument is tuned to have an unrealistically high resolution by reducing the sizes of the instrument components, such as the source, monochromator, sample and analyser. The energy FWHM of the instrument at the elastic line was found to be 0.08 meV using a simulation with the Incoherent sample component. Furthermore, the simulations do not include any sources of background or elastic scattering from the crystal. All of this aids in the comparison between the simulated datasets and the linear spin wave theory results. For the simulations, an initial number of 10⁵ rays were used for each collected data point.

The data for the simulated dispersion were obtained as constant-energy scans along four paths in q space, each of which was confined to the (h0l) plane. Each scan was performed with a step size of 0.01 meV from 0 to 8 meV for the antiferromagnet and 0 to 13.2 meV for the ferromagnet. The sample temperature was set to K for all simulations, much lower than the ordering temperature of MnF₂, $[T_{\rm N} = 67]$ K (Okazaki et al., 1964 ).

We will start by showing that the ferromagnet produces the result expected from theory. To do this, we artificially set all exchange interactions to be FM (negative), keeping their magnitude the same as those for AFM MnF₂. The exchange interactions are displayed in Table 2.

4.1. The ferromagnetic case

The ferromagnetic dispersion has been simulated along (10l). A magnetic field of T is applied to increase the size of the gap at the Brillouin zone centre. This is done to separate the spin wave creation and annihilation signals. A very small overlap between these two signals is still present, but the effect is relatively small. Gaussian line shapes were fitted to the constant-q scans to give a rough estimate of the intensity peak centres. The analytic expression for the dispersion given by equation (8) was fitted to the $[{\bf q}]$ dependence of the fitted peak positions, and the resulting parameters were compared with the original model parameters to demonstrate the consistency of our method. Since the scan is along l, only j_c, j and D were kept as free parameters in the fit (Table 3). The intensity data, fitted peak positions, and dispersions calculated from the fitted parameters and the original model parameters, respectively, can all be seen in Fig. 2 (top left). In Fig. 2 (bottom left) we show how the fitted peak positions and fitted dispersion (the dispersion calculated from the parameters extracted from the fit) deviate from the model dispersion. The deviations for both the points and the fitted dispersion are small. They generally follow each other as expected, except at q = (1, 0, 1), where mixing with the annihilation signals leads to lower energy values for the fitted peak positions.

Table 3
Ferromagnetic dispersion parameters obtained from a fit to constant-q scans, compared with the original values used in the simulation

Parameter	Value from fit [meV]	Input value to simulation [meV]
j	−0.303184 (8)	−0.304
j_c	−0.05659 (3)	−0.056
D	−±0.0332 (1)	−0.023

Figure 2
Simulated data for the ferromagnetic systems. (Top left) Simulated dispersion compared with theory. Red points show the fitted peak positions. The solid white line shows the dispersion calculated from the fitted parameters and the dashed magenta line shows the dispersion calculated from the original model parameters. (Bottom left) The difference between the fitted peak positions and the original dispersion (red points), and the difference between the fitted dispersion and the original dispersion (solid line). Ideally, they should be as close to zero as possible (dashed line). (Right) Integrated differential cross section compared with theory.

Although the shape of the fitted dispersion is seen to match theory, the extracted parameters do not exactly match those used in the simulation; the value of the anisotropy constant D deviates by around 50%. However, given that D is very small and mainly determined by the data near the gap, where the creation and annihilation signals overlap, it is not surprising that this analysis exhibits minor deviations. For the antiferromagnetic case we show how a more in-depth analysis is able to reproduce the simulation parameters with a much higher accuracy.

In addition to fitting the dispersion, we compare the simulated relative intensities with those expected from theory. The intensities are first normalized to the incoming flux and then integrated along the energy transfer, which also removes the instrument-dependent energy broadening. The values from theory are normalized to the point of the data, and the comparison is shown in Fig. 2 (right). The simulated data are seen to reproduce the theoretical values with excellent agreement.

4.2. The antiferromagnetic case

For the antiferromagnetic system, we simulate constant-q scans from the classical antiferromagnet MnF₂. The Hamiltonian input parameters were taken from Okazaki et al. (1964) and Nikotin et al. (1969 ) and are listed in Table 2(b).

The dispersion was simulated along four directions in reciprocal space, as illustrated in Fig. 3.

Figure 3
Map over scanned directions in reciprocal space for the antiferromagnetic case. Labels correspond to those in Fig. 5.

As in the ferromagnetic case, the centre of each constant-q scan intensity peak was found by fitting Gaussian line shapes to the data, as illustrated in Fig. 4. Due to the instrumental resolution function, the shape of the intensity peaks varies between different values of q. A symmetric Gaussian function is used to fit most peaks. However, near the Brillouin zone centres, at the bottom of the dispersion, the resolution function makes the peaks highly asymmetric. For these scans, an approximate resolution convoluted gap function for TAS instruments [described by Lenander et al. (2026 )] is used to fit the mode positions. In this function, it is assumed that the q resolution is coarse in two out of three directions, and at the gap, the dispersion follows a parabola, which are both valid in this setup. The function describes the asymmetric shape well; examples of the fits can be observed in Fig. 4 at low energies. Due to the extremely high energy resolution of the instrument and the assumption that the energy resolution is not correlated to the q resolution, the fitted mode positions are slightly underestimated (deviations less than 0.1 meV).

Figure 4
Examples of constant-q scans for the antiferromagnetic case (left) along (10l) and (right) along (h01). The solid lines show the fitted functions used to find the intensity peak centres and the dashed lines show the peak centres from theory.

The analytical expression for the dispersion given by equation (19) is fitted to the q dependence of the fitted peak positions. We fitted all four datasets at once, where only the interaction parameters j, j_a, j_c and D are allowed to vary. This gives us a direct comparison with the theoretical dispersion relation.

In Fig. 5 we perform this comparison by displaying the raw intensity data, the peak positions obtained from the constant-q fits, and the dispersion relations calculated from the fitted parameters and the original model parameters, respectively.

Figure 5
Simulated intensity data measured along four directions in q space. Points show fits to intensity peaks. The solid white line is the dispersion calculated from the fitted parameters and the dashed magenta line shows the dispersion calculated from the original model parameters.

The values of the fitted parameters are reported in Table 4, with the statistical uncertainties from the fits. The general agreement with the actual interaction values used as input to the component is excellent, with the largest deviation being only around 2.6% (for the small D parameter). This shows that our simulation method produces the correct dispersion and is self-consistent.

Table 4
Antiferromagnetic dispersion parameters obtained from a fit to constant-q scans, compared with the original values used in the simulations

Parameter	Value from fit [meV]	Input value to simulation [meV]
j	0.30372 (1)	0.304
j_a	0.00800 (2)	0.008
j_c	−0.05713 (2)	−0.056
D	−0.02359 (1)	−0.023

We will further validate the intensities simulated by the component used for the case of MnF₂. This is done by recording the simulated differential cross section, obtained through equation (28), at several points on the dispersion curve. In the simulations, we directly recorded the incoming flux and the outgoing intensity by placing energy-sensitive monitors just before the sample and analyser in the TAS setup. The solid angle was calculated from known geometry parameters.

We compare these values with the results of spin wave theory through equation (33), as shown in Table 5. The agreement is generally excellent. The largest relative deviation is around 12% for the very weak annihilation signal: $[\hbar\omega = -1.06]$ meV. The deviation is around 6% near the low-intensity ferromagnetic point (1, 0, 1). All other deviations are on the order of 2% or lower.

Table 5
Simulated (`sim') differential cross sections compared with results calculated from the analytic (`ana') expression at different points on the dispersion, defined by the point $[({\bf q},\hbar\omega)]$ and index of the simulated mode

Here Ψ_mon is the simulated flux at the monitor placed just before the sample, while I_det is the simulated intensity at the detector placed just before the analyser in the TAS setup.

h	k	l	ℏω [meV]	Mode	B [T]	Ψ_mon [cm⁻² s⁻¹]	I_det [s⁻¹]	(dσ/dΩ)(sim) [cm⁻²]	(dσ/dΩ)(ana) [cm⁻²]
1	0	0	1.06	0	0	$[3.509\times 10^{5}]$	$[3.750\times 10^{-1}]$	$[1.069\times 10^{-2}]$	$[1.074\times 10^{-2}]$
1	0	0	−1.06	0	0	$[2.790\times 10^{4}]$	$[5.997\times 10^{-5}]$	$[2.149\times 10^{-5}]$	$[2.453\times 10^{-5}]$
1	0	1	1.06	0	0	$[3.535\times 10^{5}]$	$[5.162\times 10^{-3}]$	$[1.460\times 10^{-4}]$	$[1.378\times 10^{-4}]$
1	0	0.5	6.73	0	0	$[5.518\times 10^{5}]$	$[6.193\times 10^{-2}]$	$[1.122\times 10^{-3}]$	$[1.097\times 10^{-3}]$
0.5	0	0.5	6.65	0	0	$[5.472\times 10^{5}]$	$[7.499\times 10^{-2}]$	$[1.370\times 10^{-3}]$	$[1.363\times 10^{-3}]$
2	0	0.3	5.48	0	0	$[5.070\times 10^{5}]$	$[1.550\times 10^{-2}]$	$[3.069\times 10^{-4}]$	$[3.030\times 10^{-4}]$
0.7	0	0.1	5.15	0	0	$[4.957\times 10^{5}]$	$[5.472\times 10^{-2}]$	$[1.104\times 10^{-3}]$	$[1.085\times 10^{-3}]$
1	0	0	1.29	0	2	$[3.609\times 10^{5}]$	$[3.882\times 10^{-1}]$	$[1.076\times 10^{-2}]$	$[1.065\times 10^{-2}]$
1	0	0	0.83	1	2	$[3.472\times 10^{5}]$	$[3.795\times 10^{-1}]$	$[1.093\times 10^{-2}]$	$[1.089\times 10^{-2}]$

The simulated intensities are also compared with a simulation of the dynamical correlation function performed using the MATLAB library SpinW (Toth & Lake, 2015). The datasets simulated along the (10l) direction were used for this. To compare the two simulated datasets, each was converted to quantities proportional to the inelastic cross section. The McStas intensities were normalized and integrated along the energy transfer, as was done for the ferromagnetic data.

The SpinW data were multiplied by the instrument-dependent factor $[k_{\rm f}/k_{\rm i}]$ , and this was then normalized to the point of the McStas data for comparison. The two datasets can be seen in Fig. 6, where we also show that our model is compatible with the SpinW simulation. An excellent agreement is found, demonstrating that the SpinWave_BCO component is able, within a constant scaling factor, to accurately reproduce results obtained by SpinW.

Figure 6
Comparison of integrated inelastic cross section between McStas, theory and SpinW. The theoretical model is shown to be compatible with the SpinW simulation. Data have been normalized to the

point.

The difference between the direct measurement of the differential cross section in Table 5 and the integrated cross section in Fig. 6 is a clear demonstration of the two distinctly different ways to integrate through reciprocal space: the differential cross section is an integration over $[E_{\rm f}]$ while keeping the scattering angle fixed, leading to the Jacobian, while Fig. 6 shows an integration over $[E_{\rm f}]$ , keeping $[{\bf q}]$ fixed, leading to data that are symmetric within a given Brillouin zone. As we demonstrate, these two equivalent methods both yield accurate results.

To further demonstrate the capabilities of the component, a simulation was performed with an applied magnetic field of T. The data are shown in Fig. 7. The two modes are split by the magnetic field, as expected. The magenta lines show the calculated spin wave mode dispersions, demonstrating that the simulated data again match well with theory. The area within the red outline was simulated using a q spacing six times smaller than the rest of the data. These data are shown integrated along the energy transfer on the right in Fig. 7. This clearly illustrates how the splitting is not always resolved. In this case this is due to the orientation of the resolution ellipsoid, leading to a focusing and a defocusing branch of the dispersion relation remeasurements simulated on the triple-axis spectrometer (Shirane et al., 2002 ), which is also seen in Fig. 5. This effect is important to consider, even for our high-resolution TAS.

Figure 7
Result of simulations made with an applied magnetic field of

T. (Left) The full $[({\bf q},\hbar\omega)]$ map showing all simulation data. (Right) A cut through these data at $[\hbar\omega = 4]$ meV and integrated along energy transfer with $[\Delta E = \pm 0.055]$ meV (red dashed outline in the left panel).

5. Altermagnetism in MnF₂

As an additional functionality, we have implemented a simple altermagnetic dispersion. Altermagnetism in MnF₂ was suggested by Šmejkal et al. (2022 ). While previous searches for the altermagnetic splitting of the magnon modes in MnF₂ were unsuccessful (Morano et al., 2025 ), a recent study revealed this splitting using polarized neutron scattering (Faure et al., 2025 ). In the latter study, the spin-wave spectrum was modelled using both altermagnetic exchange couplings and long-ranged dipolar couplings. We will use a simpler model, only including the altermagnetic couplings, to show how the SpinWave_BCO component can be used for this purpose.

The theory of exchange-driven altermagnetic spin waves in MnF₂ is described by McClarty et al. (2025 ). The lowest-order exchange interactions which lead to altermagnetic ordering are those between Mn ²⁺ ions on the same sublattice at positions characterized by the vectors $[{\hat{\bf x}}+{\hat{\bf y}}]$ and $[{\hat{\bf x}}-{\hat{\bf y}}]$ , respectively. These coupling constants will be denoted j₁₁₀ and $[j_{110}^{\,\prime}]$ , respectively. Including these in the Hamiltonian and performing the same diagonalization procedure as for the antiferromagnetic case leads to the altermagnetic dispersion relation. First, the Fourier transform of the coupling constants within each sublattice reads

$[J_{1}^{(\rm AM)}({\bf q^{\prime}}) = J_{1}({\bf q^{\prime}})+2(\,j_{110}+j^{ \,\prime}_{110})\cos\left(aq_{x}^{\prime}\right)\cos\left(bq_{y}^{\prime}\right) \eqno(34)]$

with $[J_1({\bf q^{\prime}})]$ given by equation (12). Using the new $[J_{1}^{(\rm AM)}({\bf q^{\prime}})]$ , the spin wave dispersion becomes

$[\hbar\omega^{(\rm AM)}_{{\bf q}^{\prime},a} = \Omega_{{\bf q^{\prime}}}\pm[2S(\,j_{110 }-j^{\,\prime}_{110})\sin\left(aq_{x}^{\prime}\right)\sin\left(bq_{y}^{\prime} \right)+g\mu_{\rm B}B].\eqno(35)]$

For $[j_{110} = j^{\,\prime}_{110}]$ , this reduces to the known antiferromagnetic result. For $[j_{110}\neq j^{\,\prime}_{110}]$ the two modes are split by the extra $[{\bf q^{\prime}}]$ -dependent term. This term is essentially added to the magnetic field, and just as the coherence factor is independent of the field, it is also independent of this term. Therefore, the intensity is essentially unchanged by the altermagnetic splitting. For this reason, we have implemented the altermagnetic splitting by simply modifying the dispersion relation according to equation (35). A simulation of the altermagnetic dispersion has been performed on the same TAS setup as used before.

To show the splitting, the values of j₁₁₀ and $[j_{110}^{\,\prime}]$ and the difference between them are made very large: $[j_{110} = -0.25 j]$ and $[j_{110}^{\,\prime} = -0.5 j]$ . The simulated spectrum can be seen in Fig. 8, where the theoretical dispersion relation has been plotted for comparison. The altermagnetic splitting is clearly seen in the data.

Figure 8
Simulated altermagnetic dispersion along (h,h+1,0). The dashed magenta lines show the theoretical dispersion relations for the split spin wave modes.

6. Discussion and conclusion

We have demonstrated the implementation of accurate magnetic dynamics in McStas with a new sample component SpinWave_BCO, which simulates inelastic neutron scattering from spin waves in a ferromagnet or a two-sublattice antiferromagnet in a body-centred orthorhombic crystal structure. In addition, we have implemented and tested an altermagnetic splitting of the antiferromagnetic spin wave dispersion.

The simulated dispersion relation and differential cross section are found from linear spin wave theory. Simulated data have been compared directly with theory and with the program SpinW to verify the McStas simulations. This was exemplified by simulating the antiferromagnetic spin wave spectrum of MnF₂. However, since the magnetic form factor is not taken into account, the simulated intensities are reliable only at small values of |q|.

We find that the simulated dispersions and absolute intensities agree extremely well with linear spin wave theory, both for magnon creation and for annihilation. This has been tested for many different values of |q|.

While our present work is an important step towards the description of magnetic dynamics in McStas, much development work lies ahead in this direction. The feature of greatest relevance will surely be to expand the description of lattice geometries and interactions to be able to model all crystal classes and also accommodate spin waves from n-sublattice systems and spiral order, in analogy with SpinW (Toth & Lake, 2015). To promote the calculation of realistic intensities, it could be of relevance to include a library of the magnetic structure factors, $[F({\bf q})]$ , for the most common magnetic ions. It would also be relevant to evaluate other possible sampling methods, such as a pre-compiled grid of solutions to the kinematic constraint (26). A very ambitious, but also relevant, update would be to include a full description of neutron polarization, which is generally accommodated in McStas (Knudsen et al., 2014 ) but the effect of which is not included in the SpinWave_BCO component.

Acknowledgements

We thank Rasmus Toft-Petersen for expressing the need for this spin wave component to be developed and Henrik M. Rønnow for suggesting to include the altermagnetic dispersion. Furthermore, we thank Amalie F. Davidsen and Frida B. Nielsen for useful discussions on inelastic scattering in McStas. We also thank Jakob Lass for supplying the SpinW simulations of MnF₂ used to validate the McStas simulations. Finally, we thank Peter Willendrup and Mads Bertelsen for including the SpinWave_BCO component in the McStas component library.

Conflict of interest

The authors declare no conflicts of interest.

Data availability

All simulation files, simulated data and scripts for data analysis can be obtained from the authors upon reasonable request. The SpinWave_BCO component is included in the McStas component library, available from the homepage https://www.mcstas.org/.

Funding information

This project was supported by the Danish Committee for Research Infrastructure (NUFI) through funding to the `ESS-Lighthouse' Q-MAT.

References

Andersen, K. H., Argyriou, D. N., Jackson, A. J., Houston, J., Henry, P. F., Deen, P. P., Toft-Petersen, R., Beran, P., Strobl, M., Arnold, T., Wacklin-Knecht, H., Tsapatsaris, N., Oksanen, E., Woracek, R., Schweika, W., Mannix, D., Hiess, A., Kennedy, S., Kirstein, O., Petersson Årsköld, S., Taylor, J., Hagen, M. E., Laszlo, G., Kanaki, K., Piscitelli, F., Khaplanov, A., Stefanescu, I., Kittelmann, Th., Pfeiffer, D., Hall-Wilton, R., Lopez, C. I., Aprigliano, G., Whitelegg, L., Moreira, F. Y., Olsson, M., Bordallo, H. N., Martín-Rodríguez, D., Schneider, H., Sharp, M., Hartl, M., Nagy, G., Ansell, S., Pullen, S., Vickery, A., Fedrigo, A., Mezei, F., Arai, M., Heenan, R. K., Halcrow, W., Turner, D., Raspino, D., Orszulik, A., Cooper, J., Webb, N., Galsworthy, P., Nightingale, J., Langridge, S., Elmer, J., Frielinghaus, H., Hanslik, R., Gussen, A., Jaksch, S., Engels, R., Kozielewski, T., Butterweck, S., Feygenson, M., Harbott, P., Poqué, A., Schwaab, A., Lieutenant, K., Violini, N., Voigt, J., Brückel, T., Koenen, M., Kämmerling, H., Babcock, E., Salhi, Z., Wischnewski, A., Heynen, A., Désert, S., Jestin, J., Porcher, F., Fabrèges, X., Fabrèges, G., Annighöfer, B., Klimko, S., Dupont, Th., Robillard, Th., Goukassov, A., Longeville, S., Alba-Simionesco, Ch., Bourges, P., Guyon Le Bouffy, J., Lavie, P., Rodrigues, S., Calzada, E., Lerche, M., Schillinger, B., Schmakat, P., Schulz, M., Seifert, M., Lohstroh, W., Petry, W., Neuhaus, J., Loaiza, L., Tartaglione, A., Glavic, A., Schütz, S., Stahn, J., Lehmann, E., Morgano, M., Schefer, J., Filges, U., Klauser, Ch., Niedermayer, Ch., Fenske, J., Nowak, G., Rouijaa, M., Siemers, D. J., Kiehn, R., Müller, M., Carlsen, H., Udby, L., Lefmann, K., Birk, J. O., Holm-Dahlin, S., Bertelsen, M., Hansen, U. B., Olsen, M. A., Christensen, M., Iversen, K., Christensen, N. B., Rønnow, H. M., Freeman, P. G., Hauback, B. C., Kolevatov, R., Llamas-Jansa, I., Orecchini, A., Sacchetti, F., Petrillo, C., Paciaroni, A., Tozzi, P., Zanatta, M., Luna, P., Herranz, I., del Moral, O. G., Huerta, M., Magán, M., Mosconi, M., Abad, E., Aguilar, J., Stepanyan, S., Bakedano, G., Vivanco, R., Bustinduy, I., Sordo, F., Martínez, J. L., Lechner, R. E., Villacorta, F. J., Šaroun, J., Lukáš, P., Markó, M., Zanetti, M., Bellissima, S., del Rosso, L., Masi, F., Bovo, C., Chowdhury, M., De Bonis, A., Di Fresco, L., Scatigno, C., Parker, S. F., Fernandez-Alonso, F., Colognesi, D., Senesi, R., Andreani, C., Gorini, G., Scionti, G. & Schreyer, A. (2020). Nucl. Instrum. Methods Phys. Res. A 957, 163402. Web of Science CrossRef Google Scholar
Boothroyd, A. T. (2020). Principles of Neutron Scattering from Condensed Matter. Oxford, New York: Oxford University Press. Google Scholar
Christensen, D. L. et al. (2026). J. Appl. Cryst. In preparation. Google Scholar
Davidsen, A. F., Nielsen, F. B., Krighaar, K. M. L., Wang, X., Avery, J. & Lefmann, K. (2026). J. Appl. Cryst. In preparation. Google Scholar
Faure, Q., Bounoua, D., BaléDent, V., Gukasov, A., Garlea, V. O., Ribeiro, A., Rau, J. G., Petit, S. & McClarty, P. (2025). arXiv, 2509.07087, https://arxiv.org/abs/2509.07087. Google Scholar
Knudsen, E. B., Tranum-Rømer, A., Willendrup, P. K., Christiansen, P. & Lefmann, K. (2014). J. Neutron Res. 17, 27–34. CrossRef CAS Google Scholar
Lefmann, K. & Nielsen, K. (1999). Neutron News 10(3), 20–23. CrossRef Google Scholar
Lefmann, K., Willendrup, P. K., Udby, L., Lebech, B., Mortensen, K., Birk, J. O., Klenø, K., Knudsen, E., Christiansen, P., Saroun, J., Kulda, J., Filges, U., Konnecke, M., Tregenna-Piggott, P., Peters, J., Lieutenant, K., Zsigmond, G., Bentley, P. & Farhi, E. (2008). J. Neutron Res. 16, 97–111. CrossRef CAS Google Scholar
Lenander, E. Y., Schack, S. B., Lefmann, K. & Rønnow, H. M. (2026). J. Appl. Cryst. 59, 783–798. CrossRef IUCr Journals Google Scholar
Lin, J. Y. Y., Islam, F., Sala, G., Lumsden, I., Smith, H., Doucet, M., Stone, M. B., Abernathy, D. L., Ehlers, G., Ankner, J. F. & Granroth, G. E. (2019). J. Phys. Commun. 3, 085005. CrossRef Google Scholar
Lin, J. Y. Y., Smith, H. L., Granroth, G. E., Abernathy, D. L., Lumsden, M. D., Winn, B., Aczel, A. A., Aivazis, M. & Fultz, B. (2016). Nucl. Instrum. Methods Phys. Res. A 810, 86–99. Web of Science CrossRef CAS Google Scholar
Marshall, W. & Lovesey, S. W. (1971). Theory of Thermal Neutron Scattering. Oxford University Press. Google Scholar
McClarty, P. A., Gukasov, A. & Rau, J. G. (2025). Phys. Rev. B 111, L060405. Web of Science CrossRef Google Scholar
Morano, W. C., Maesen, Z., Nikitin, S. E., Lass, J., Mazzone, D. G. & Zaharko, O. (2025). Phys. Rev. Lett. 134, 226702. Web of Science CrossRef PubMed Google Scholar
Nikotin, O., Lindgård, P. A. & Dietrich, O. W. (1969). J. Phys. C Solid State Phys. 2, 1168–1173. CrossRef CAS Web of Science Google Scholar
Okazaki, A., Turberfield, K. & Stevenson, R. (1964). Phys. Lett. 8, 9–11. CrossRef CAS Web of Science Google Scholar
Shirane, G., Shapiro, S. M. & Tranquada, J. M. (2002). Neutron Scattering with a Triple-Axis Spectrometer. Cambridge University Press. Google Scholar
Šmejkal, L., Sinova, J. & Jungwirth, T. (2022). Phys. Rev. X 12, 031042. Google Scholar
Squires, G. L. (1977). Introduction to the Theory of Thermal Neutron Scattering. Cambridge University Press. Google Scholar
Toth, S. & Lake, B. (2015). J. Phys. Condens. Matter 27, 166002. Web of Science CrossRef PubMed Google Scholar
Willendrup, P., Filges, U., Keller, L., Farhi, E. & Lefmann, K. (2006). Physica B 385–386, 1032–1034. Web of Science CrossRef CAS Google Scholar
Willendrup, P. K. & Lefmann, K. (2020). J. Neutron Res. 22, 1–16. Web of Science CrossRef Google Scholar
Willendrup, P. K. & Lefmann, K. (2021). J. Neutron Res. 23, 7–27. Web of Science CrossRef Google Scholar