# Magnetic Neutron Scattering from Spherical Nanoparticles with Néel Surface Anisotropy: Analytical Treatment

Michael P. Adams, Andreas Michels, Hamid Kachkachi

► To cite this version:

Michael P. Adams, Andreas Michels, Hamid Kachkachi. Magnetic Neutron Scattering from Spherical Nanoparticles with Néel Surface Anisotropy: Analytical Treatment. Journal of Applied Crystallography, 2022, 55 (6), pp.1475-1487. 10.1107/S1600576722008925. hal-04751007

HAL Id: hal-04751007

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

Submitted on 25 Oct 2024

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>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

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

JOURNAL OF
APPLIED
CRYSTALLOGRAPHY

ISSN 1600-5767

# Magnetic neutron scattering from spherical nanoparticles with Néel surface anisotropy: analytical treatment

Michael P. Adams,<sup>a,*</sup> Andreas Michels<sup>a</sup> and Hamid Kachkachi<sup>b</sup>

<sup>a</sup>Department of Physics and Materials Science, University of Luxembourg, 162A avenue de la Faiencerie, L-1511 Luxembourg, Grand Duchy of Luxembourg, and <sup>b</sup>Laboratoire PROMES CNRS UPR8521, Université de Perpignan via Domitia, Rambla de la Thermodynamique, Tecnosud, F-66100 Perpignan, France. *Correspondence e-mail: michael.adams@uni.lu

The magnetization profile and the related magnetic small-angle neutron scattering cross section of a single spherical nanoparticle with Néel surface anisotropy are analytically investigated. A Hamiltonian is employed that comprises the isotropic exchange interaction, an external magnetic field, a uniaxial magnetocrystalline anisotropy in the core of the particle and the Néel anisotropy at the surface. Using a perturbation approach, the determination of the magnetization profile can be reduced to a Helmholtz equation with Neumann boundary condition, whose solution is represented by an infinite series in terms of spherical harmonics and spherical Bessel functions. From the resulting infinite series expansion, the Fourier transform, which is algebraically related to the magnetic small-angle neutron scattering cross section, is analytically calculated. The approximate analytical solution for the spin structure is compared with the numerical solution using the Landau–Lifshitz equation, which accounts for the full nonlinearity of the problem. The signature of the Néel surface anisotropy can be identified in the magnetic neutron scattering observables, but its effect is relatively small, even for large values of the surface anisotropy constant.

## 1. Introduction

Magnetic small-angle neutron scattering (SANS) is a powerful technique for investigating spin structures on the mesoscopic length scale (~1–100 nm) and inside the volume of magnetic materials (Mühlbauer et al., 2019; Michels, 2021). Recent SANS studies of magnetic nanoparticles, in particular employing spin-polarized neutrons, demonstrate that their spin textures are highly complex and exhibit a variety of nonuniform, canted or core–shell-type configurations [see e.g. Disch et al. (2012), Krycka et al. (2014), Hasz et al. (2014), Günther et al. (2014), Maurer et al. (2014), Dennis et al. (2015), Grutter et al. (2017), Oberdick et al. (2018), Ijiri et al. (2019), Bender et al. (2019), Bersweiler et al. (2019), Zákutná et al. (2020), Honecker et al. (2022) and references therein]. Surface anisotropy, vacancies or the presence of antiphase boundaries are generally considered to be at the origin of spin disorder in nanoparticles (Berger et al., 2008; Wetterskog et al., 2013; Nedelkoski et al., 2017; Köhler et al., 2021a; Batlle et al., 2022). Magnetic SANS data analysis largely relies on structural form-factor models for the cross section, borrowed from nuclear SANS, which do not properly account for the existing spin inhomogeneity inside magnetic nanoparticles or nanomagnets (NMs).

Received 16 May 2022
Accepted 6 September 2022

Edited by G. J. McIntyre, Australian Nuclear Science and Technology Organisation, Lucas Heights, Australia

**Keywords:** magnetic neutron scattering; small-angle neutron scattering; magnetic nanoparticles; surface anisotropy; micromagnetics.

**Supporting information:** this article has supporting information at journals.iucr.org/j

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

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

OPEN

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

ACCESS

Published under a CC BY 4.0 licence
Progress in magnetic SANS theory (Honecker & Michels, 2013; Michels et al., 2014; Mettus & Michels, 2015; Erokhin et al., 2015; Metlov & Michels, 2015, 2016; Michels et al., 2016, 2019; Mistonov et al., 2019; Zaporozhets et al., 2022) strongly suggests that, for the analysis of experimental magnetic SANS data, the spatial nanometre-scale variation of the orientation and magnitude of the magnetization vector field must be taken into account and macrospin-based models – assuming a uniform magnetization – are not adequate. The starting point for a proper analysis of the scattering problem is a micro-magnetic continuum expression for the magnetic energy of the system. In the static case, this then leads to Brown's equations (Brown, 1963), a set of nonlinear partial differential equations for the magnetization along with complex boundary conditions on the surface of the magnet. From these equations the Fourier image and the magnetic SANS cross section may be obtained.

In this paper, we present an analytical treatment of the magnetic SANS cross section of a spherical NM with Néel surface anisotropy (Néel, 1954). This particular form of anisotropy arises because in an NM a significant fraction of atoms belong to the surface (with no neighbours on one side), and their magnetic properties such as exchange and anisotropy can be strongly modified relative to the bulk atoms.

The manuscript is organized as follows. In Section 2, we calculate the real-space spin structure of a spherical NM using classical micromagnetic theory within the second-order perturbation approach. In Section 3, we compute the three-dimensional Fourier transform of the real-space spin structure, which directly yields the magnetic neutron scattering cross section and the pair-distance distribution function. The analytical results are benchmarked by comparing them with numerical finite difference simulations using the Landau-Lifshitz equation of motion. Finally, Section 4 summarizes the main findings of this study.

We also make reference to our accompanying numerical study (Adams et al., 2022) where, in contrast to the present analytical work, the full nonlinearity of the problem is considered.

## 2. Micromagnetic theory

In the static micromagnetic approach (Brown, 1963), the magnetic configuration of a system is described by the continuous magnetization vector field $M(r)$, which has a constant magnitude $||M(r)|| = M_0$. The saturation magnetization $M_0$ is only a function of temperature. The normalized magnetization vector field is then defined as

$$ \mathbf{m}(\mathbf{r}) = \mathbf{M}(\mathbf{r})/M_0 = [m_x(\mathbf{r}), m_y(\mathbf{r}), m_z(\mathbf{r})], \quad (1) $$

where $\mathbf{r}$ denotes the position vector. Our Hamiltonian for the NM includes the isotropic exchange interaction, the Zeeman energy, a uniaxial magnetic anisotropy for spins in the core and Néel surface anisotropy for those on the surface. In the continuum approach, it reads

$$ \begin{aligned} \mathcal{H} = & -A \sum_{\alpha \in \{x,y,z\}} \int_V m_{\alpha} \Delta m_{\alpha} \, d^3 r - M_0 \mathbf{B}_0 \int_V \mathbf{m} \, d^3 r \\ & - K_c \int_V (\mathbf{m} \cdot \mathbf{e}_A)^2 \, d^3 r + A \sum_{\alpha \in \{x,y,z\}} \oint_{\partial V} m_{\alpha} \nabla m_{\alpha} \cdot \mathbf{n} \, d^2 r \\ & - \frac{K_s}{2} \sum_{\alpha \in \{x,y,z\}} \oint_{\partial V} |n_{\alpha}| m_{\alpha}^2 \, d^2 r, \end{aligned} \quad (2) $$

where $A$ is the exchange stiffness constant, $\nabla$ is the del operator, $\Delta$ is the Laplace operator, $\mathbf{B}_0$ is a constant applied magnetic field, $K_c > 0$ denotes the uniaxial core anisotropy constant, $\mathbf{e}_A$ is a unit vector specifying the arbitrary core anisotropy axis and $K_s > 0$ is the Néel surface anisotropy constant (Néel, 1954).

$$ \mathbf{n} = [\sin \theta \cos \phi, \sin \theta \sin \phi, \cos \theta] \quad (3) $$

is the surface normal to the boundary of the NM, where $\theta$ and $\phi$ are the usual spherical angles (Garanin & Kachkachi, 2003; Kachkachi, 2007). In (2), the two surface integrals take into account the boundary conditions for the magnetization on the surface ($\partial V$) of the NM of volume $V$, which result from the exchange interaction and the Néel term. The magnetodipolar energy has been ignored in the calculations because of its mathematical complexity and since it is expected to be of minor relevance for smaller-sized NMs [see the recent atomistic simulations by Köhler et al. (2021b)].

For small deviations from the homogeneous magnetization state, a perturbation approach is applicable. Let $\mathbf{m}_0$ be the principal unit vector (average direction) associated with $\mathbf{m}(\mathbf{r})$ and let the vector function $\psi(\mathbf{r}) \perp \mathbf{m}_0$ describe the spin misalignment. One can then write

$$ \mathbf{m}(\mathbf{r}) = \mathbf{m}_0 [1 - \|\psi(\mathbf{r})\|^2]^{1/2} + \psi(\mathbf{r}), \quad \|\mathbf{m}(\mathbf{r})\| = 1. \quad (4) $$

Assuming that $\psi_x, \psi_y, \psi_z \ll 1$, the following second-order Maclaurin expansion in $\psi$ is used to find an approximate closed-form solution for $\mathbf{m}(\mathbf{r})$:

$$ \mathbf{m}(\mathbf{r}) \approx \mathbf{m}_0 + \psi(\mathbf{r}) - \frac{1}{2} \|\psi(\mathbf{r})\|^2 \mathbf{m}_0, \quad (5) $$

where $\mathbf{m}_0$ is taken as a known constant vector in subsequent calculations. We choose the orthonormal vector base (Garanin & Kachkachi, 2009),

$$ \left\{ \begin{array}{l} \mathbf{g}_0 = \mathbf{m}_0, \\ \mathbf{g}_1 = \frac{\mathbf{m}_0 \times \mathbf{e}_A}{\|\mathbf{m}_0 \times \mathbf{e}_A\|}, \\ \mathbf{g}_2 = \frac{(\mathbf{m}_0 \cdot \mathbf{e}_A)\mathbf{m}_0 - \mathbf{e}_A}{\|(\mathbf{m}_0 \cdot \mathbf{e}_A)\mathbf{m}_0 - \mathbf{e}_A\|}, \end{array} \right. \quad (6) $$

and the parametrization

$$ \psi(\mathbf{r}) = \psi_1(\mathbf{r})\mathbf{g}_1 + \psi_2(\mathbf{r})\mathbf{g}_2, \quad (7) $$

and introduce the dimensionless coordinates $\xi = \mathbf{r}/R$ (with $\xi = \|\xi\| = r/R$), where $\mathbf{r}$ is the position vector,

$$ \mathbf{r} = [r \sin \theta \cos \phi, r \sin \theta \sin \phi, r \cos \theta], \quad (8) $$

and $R$ denotes the radius of the NM. The minimization of the Hamiltonian (2) then leads to the well known Helmholtz
equation with Neumann boundary conditions on the unit sphere (Kachkachi, 2007; Garanin & Kachkachi, 2003),

$$[\Delta_{\xi} - \kappa_{\beta}^{2}] \psi_{\beta} = 0, \quad \beta \in \{1, 2\}, \qquad (9)$$

$$\left. \frac{d\psi_{\beta}}{d\xi} \right|_{\xi=1} = \sum_{\alpha \in \{x, y, z\}} \chi_{\alpha}^{\beta} |n_{\alpha}|, \qquad (10)$$

where the constants are defined as

$$\kappa_{1}^{2} = \mathbf{m}_{0} \cdot \mathbf{b}_{0} + 2k_{c}(\mathbf{m}_{0} \cdot \mathbf{e}_{\Lambda})^{2}, \qquad (11)$$

$$\kappa_{2}^{2} = \mathbf{m}_{0} \cdot \mathbf{b}_{0} + 2k_{c}[2(\mathbf{m}_{0} \cdot \mathbf{e}_{\Lambda})^{2} - 1], \qquad (12)$$

$$\chi_{\alpha}^{\beta} = k_{s}(\mathbf{m}_{0} \cdot \mathbf{e}_{\alpha}) (\mathbf{g}_{\beta} \cdot \mathbf{e}_{\alpha}), \qquad (13)$$

with the dimensionless quantities

$$k_c = \frac{R^2 K_c}{2A}, \quad k_s = \frac{R K_s}{2A}, \quad \mathbf{b}_0 = \frac{R^2 M_0}{2A} \mathbf{B}_0. \qquad (14)$$

The **e**_α (with α = x, y, z) in (13) denote the unit vectors of the Cartesian laboratory coordinate frame (in which **n** and **r** are defined). We emphasize that there are only two independent differential equations for **ψ**, which is a consequence of the constraint ||**m**(**r**)|| = 1.

In our graphical representations, we will frequently use the following values: $k_c = 0.1$ and $k_s = 3.0$, which (using $R = 5$ nm and $A = 10^{-11}$ J m⁻¹) correspond to $K_c = 80$ kJ m⁻³ and $K_s = 12$ mJ m⁻² (Gradmann, 1986; O’Handley, 2000; Battlle et al., 2022). For $M_0 = 1.7 \times 10^6$ A m⁻¹, the relation between $b_0$ (dimensionless) and the external field is $B_0 = (8/17)b_0 \times 1$ T.

The fundamental solution of the homogeneous Helmholtz equation (9) is well known (Weber & Arfken, 2003; Riley et al., 2006). Its non-singular part can be expressed in spherical coordinates as an infinite series in terms of spherical harmonics $Y_{\ell m}(\theta, \phi)$ and spherical Bessel functions of the first kind $j_n(i\kappa\beta\xi)$,

$$\psi_{\beta} = \sum_{\ell=0}^{\infty} \sum_{m=-\ell}^{\ell} c_{\ell m}^{\beta} j_{\ell}(i\kappa\beta\xi) Y_{\ell m}(\theta, \phi). \qquad (15)$$

The imaginary number ‘i’ in the argument of the spherical Bessel function is due to the negative sign in the Helmholtz equation (9). The expansion coefficients $c_{\ell m}^{\beta}$ are obtained from the Neumann boundary condition (10) using the method of least squares (see Appendix A). From Appendix A it is seen that the zero-order term with $\ell = 0$ vanishes. This physically makes sense, since the spin misalignment in our model is caused by the Néel surface anisotropy and thus, for symmetry reasons, there is no misalignment at the centre of the NM, i.e. $\psi_{\beta}(\xi = 0, \theta, \phi) \equiv 0$. By contrast, the largest spin misalignment is found at the boundary of the NM, i.e. $\xi = 1$. Further, we find that the coefficients $c_{\ell m}^{\beta}$ vanish in the case of odd $\ell$ and $m$, while they are real valued and even with respect to the index $m$, i.e. $c_{\ell m}^{\beta} = c_{\ell, -m}^{\beta}$. Taking these properties into account, one can conveniently express the solution in terms of the associated Legendre polynomials $P_{\ell}^{m}(\cos\theta)$ with $\ell = 2\nu$ and $m = 2\mu$ (note that we use the convention that $Y_{\ell m}(\theta, \phi) = N_{\ell m}P_{\ell}^{m}(\cos\theta)\exp(im\phi)$ [p. 378 (14.30.1) of Olver et al. (2010)]),

$$\psi_{\beta} = \sum_{\nu=1}^{\infty} \sum_{\mu=0}^{\nu} a_{\nu\mu}^{\beta} \Upsilon_{\nu}(\kappa_{\beta}\xi) P_{2\nu}^{2\mu}(\cos\theta) \cos(2\mu\phi), \qquad (16)$$

where we define [compare pp. 624–626 of Weber & Arfken (2003)]

$$\Upsilon_{\nu}(\tau) = j_{2\nu}(i\tau) = \frac{\pi^{1/2}}{2} \sum_{s=0}^{\infty} \frac{(-1)^{\nu}(\tau/2)^{2(s+\nu)}}{s!\Gamma(2\nu + s + 3/2)}, \qquad (17)$$

and the expansion coefficients are given by

$$a_{\nu\mu}^{\beta} = \frac{2k_s N_{2\nu,2\mu} \mathbf{g}_{\beta} \cdot \text{diag}[I_{2\nu,2\mu}^{x}, I_{2\nu,2\mu}^{y}, I_{2\nu,2\mu}^{z}] \cdot \mathbf{m}_{0}}{1 + \delta_{\mu,0} \kappa_{\beta} \Upsilon'_{\nu}(\kappa_{\beta})}, \qquad (18)$$

with

$$I_{\ell m}^{\alpha} = \int_{0}^{2\pi} \int_{0}^{\pi} Y_{\ell, m}^{*}(\theta, \phi) |n_{\alpha}| \sin\theta \, d\theta \, d\phi \qquad (19)$$

and

$$N_{\ell m} = \left[ \frac{2\ell + 1}{4\pi} \frac{(\ell - m)!}{(\ell + m)!} \right]^{1/2} \qquad (20)$$

In (18), $\delta_{\mu,0}$ is the Kronecker delta function, diag[...] denotes a $3 \times 3$ diagonal matrix and $\Upsilon_{\nu}(\tau)$ is the first-order derivative of (17) with respect to $\tau$. For some small values of $\ell$ and $m$, the exact solutions of the integrals $I_{\ell m}^{\alpha}$ are listed in Table 1 in Appendix B.

From (18) it is seen that the functions $\psi_{\beta}$ depend linearly on $k_s$, such that for $k_s = 0$ the magnetization of the NM is homogeneous (as expected). Since we assume that $\psi_x, \psi_y, \psi_z \ll 1$, it is clear that the validity of our solution is restricted to a finite range $0 \le k_s \le k_{s, \max}$. Taking only the terms with $\nu = 1$ into account (corresponding to $\ell = 2$), the remaining (second-order) expression reads

$$\psi_{\beta}(\xi, \theta, \phi) \approx -\frac{15k_s}{32} \frac{\Upsilon_1(\kappa_{\beta}\xi)}{\kappa_{\beta} \Upsilon'_1(\kappa_{\beta})} \mathbf{g}_{\beta} \cdot \mathcal{V}(\theta, \phi) \cdot \mathbf{m}_0, \qquad (21)$$

where

$$\mathcal{V}(\theta, \phi) = \text{diag} \left[ \begin{array}{c} \cos^2 \theta - 1/3 - \sin^2 \theta \cos(2\phi) \\ \cos^2 \theta - 1/3 + \sin^2 \theta \cos(2\phi) \\ -2(\cos^2 \theta - 1/3) \end{array} \right]. \qquad (22)$$

A reasonable approximation for small $\kappa_{\beta}$ in (21) is obtained by taking into account the first two terms in the infinite series (17) for $\Upsilon_{\nu}(\tau)$. This results in the following expression [compare (21)]:

$$\frac{\Upsilon_1(\kappa_\beta \xi)}{\kappa_\beta \Upsilon'_1(\kappa_\beta)} \approx \frac{1}{4} \frac{\kappa_\beta^2 \xi^4 + 14 \xi^2}{\kappa_\beta^2 + 7}. \qquad (23)$$

In the limit $\kappa_{\beta} \to 0$ (small $\mathbf{B}_0$ and small $K_c$), this expression reduces to a quadratic function in $\xi$,

$$\lim_{\kappa_{\beta} \to 0} \left\{ \frac{\Upsilon_1(\kappa_{\beta}\xi)}{\kappa_{\beta}\Upsilon'_1(\kappa_{\beta})} \right\} = \frac{\xi^2}{2}. \qquad (24)$$

The case of an infinite applied magnetic field $\mathbf{B}_0$, or of a strong uniaxial core anisotropy [compare (11) and (12)], corresponds to the limit
<figure><img src="image_7.png" /><figcaption>Figure 1</figcaption></figure>

The normalized effective energy potential of the Néel surface anisotropy as a function of the Cartesian components of the average magnetization vector **m**₀ = [sin βcos α, sin βsin α, cos β], computed via numerical integration of the surface contribution in (2) and using the second-order approximation (21). Parameters are **e**<sub>A</sub> = [0, 0, 1], **b**₀ = **0**, k<sub>c</sub> = 0.1 and k<sub>s</sub> = 3.0. The minima of the Néel surface contribution are in this case along the cubic space diagonals **m**₀ = [±1, ±1, ±1]/(3<sup>1/2</sup>), while the maxima correspond to the Cartesian axes ± **e**<sub>x</sub>, ± **e**<sub>y</sub>, ± **e**<sub>z</sub>. The effective energy potential has cubic symmetry and is approximately proportional to a function of the type ≈ **m**<sub>0,x</sub><sup>4</sup> + **m**<sub>0,y</sub><sup>4</sup> + **m**<sub>0,z</sub><sup>4</sup> [see also Garanin & Kachkachi (2003)].

$$ \lim_{\kappa_{\beta} \to \infty} \left\{ \frac{\Upsilon_1(\kappa_{\beta}\xi)}{\kappa_{\beta}\Upsilon_1(\kappa_{\beta})} \right\} = 0, \qquad (25) $$

which recovers the expected result of zero spin misalignment. Note that the limit $κ_β → ∞$ is only obtained using all terms of the infinite series (17).

Of particular interest is the behaviour of $ψ_β$ as a function of the radius R of the NM. Inspecting the Hamiltonian (2), it becomes clear that the surface anisotropy energy scales as R<sup>2</sup>, while the uniaxial core anisotropy energy scales as R<sup>3</sup>. Since the core and surface anisotropies act in opposite ways (trying to make the spin structure more homogeneous and more inhomogeneous, respectively), we see that an increasing radius R corresponds to a decreasing $ψ_β$. This behaviour reflects the NM's surface-area-to-volume ratio. With (21) it is not possible to make any prediction in this regard, because until this point we have not included the principal unit vector **m**₀ in the minimization of the Hamiltonian. Generally, **m**₀ is a function of k<sub>s</sub>, k<sub>c</sub>, **b**₀ and **e**<sub>A</sub>.

In the special case when the uniaxial anisotropy axis and the applied magnetic field are both directed parallel to the z axis (e<sub>A</sub> = [0, 0, 1] and **b**₀ = [0, 0, b₀]), the principal unit magnetization vector may be written as

$$ \mathbf{m}_0 = [(1/2^{1/2}) \sin \beta, (1/2^{1/2}) \sin \beta, \cos \beta], \quad (26) $$

where $β ∈ [0, arccos(1/3^{1/2})]$. This choice is justified by the effective cubic symmetry of the Néel anisotropy as shown in

<figure><img src="image_8.png" /><figcaption>Figure 2</figcaption></figure>

A comparison between the numerical solution using the Landau–Lifshitz equation (upper row) and the second-order analytical solution (21) for $||ψ(ξ)|| = [ψ_1^2(ξ) + ψ_2^2(ξ)]^{1/2}$ (lower row). (a) and (e) show $||ψ||$ on the boundary surface ($ξ = 1$), while (b)–(d) and (f)–(h) display selected planar cuts in ($ξ_x, ξ_y, ξ_z$) space. The following parameters are used: **e**<sub>A</sub> = **e**<sub>z</sub>, **b**₀ = [0.4, 0, 0.4] (B₀ ≈ 266 mT), k<sub>c</sub> = 0.1, k<sub>s</sub> = 3.0 and **m**₀ = [sin βcos α, sin βsin α, cos β], where α = 0° and β = 40°.
Fig. 1. This result was already predicted by Garanin & Kachkachi (2003). The solutions for $\psi_{1,2}(\xi, \theta, \phi)$ [using the particular $\mathbf{m}_0$ (26)] then read

$$ \psi_1 \approx \frac{15k_s}{32} \frac{\Upsilon_1(\kappa_1\xi)}{\kappa_1 \Upsilon'_1(\kappa_1)} \sin^2\theta \cos(2\phi) \sin\beta, \quad (27) $$

$$ \psi_2 \approx \frac{15k_s}{32} \frac{\Upsilon_1(\kappa_2\xi)}{\kappa_2 \Upsilon'_1(\kappa_2)} (1 - 3\cos^2\theta) \sin\beta \cos\beta. \quad (28) $$

In Fig. 2, the analytical solution (21) (lower row) is compared with the numerical solution based on the Landau–Lifshitz equation $\dot{\mathbf{m}} = -\gamma\mathbf{m} \times \mathbf{B}_{\text{eff}} - \alpha\mathbf{m} \times (\mathbf{m} \times \mathbf{B}_{\text{eff}})$ (upper row) (Bertotti, 1998), where $\gamma$ is the gyromagnetic ratio, $\alpha$ is the damping constant and the dot denotes the first-order time derivative [see our numerical study in the accompanying paper (Adams et al., 2022) for further details]. Shown is the vector norm of the $\psi(\xi)$ function scaled to its maximum value. From Fig. 2 it is seen that our analytical approximation is in qualitative agreement with the results from the numerical simulation. The corresponding real-space spin structure $\mathbf{m}(\xi)$ is displayed in Fig. 3, where the surface spin disorder becomes clearly visible.

It is also instructive to compare our solution (21) with that obtained using the Green's function approach (Garanin & Kachkachi, 2003; Kachkachi, 2007). In particular, for $\xi$ located close to the surface, where the maximum spin misalignment with respect to $\mathbf{m}_0$ occurs, the Green's function method yields the following approximate expression:

$$ \psi_{\beta}(\xi) \approx -\frac{15k_s}{32} \left(1 - \frac{\kappa_{\beta}^2}{14}\right) \xi^2 \mathbf{g}_{\beta} \cdot \mathcal{V}(\xi) \cdot \mathbf{m}_0, \quad (29) $$

<figure><img src="image_9.png" /><figcaption>Figure 3 The real-space spin structure in the $\xi_x\xi_z$ plane computed using (4) and (21). Parameters are the same as in Fig. 2. The external field $B_0 \approx 266$ mT is applied in the $\xi_x\xi_z$ plane and inclined by an angle of $\beta = 40^\circ$ relative to the $\xi_z$ axis [compare with Garanin & Kachkachi (2003)].</figcaption></figure>

$$ \mathcal{V}(\xi) = -\operatorname{diag} \begin{bmatrix} \xi_x^2/\xi^2 - 1/3 \\ \xi_y^2/\xi^2 - 1/3 \\ \xi_z^2/\xi^2 - 1/3 \end{bmatrix}. \quad (30) $$

This expression is also found when (21) is expanded in $\kappa_\beta$ at the surface of the NM ($\xi = 1$).

While the infinite series approach using spherical harmonics and spherical Bessel functions yields an exact solution of the Helmholtz equation, the Green's function approach provides an approximate explicit expression of $\psi_\beta$ in terms of the coefficients $\kappa_\beta$. Indeed, as was shown by Kachkachi (2007), in the presence of core anisotropy Green's function as the kernel of the Helmholtz equation is only obtained as a perturbative series in $\kappa_\beta$. As such, (29) is restricted to small values of $\kappa_\beta$, i.e. assuming that the core anisotropy and applied magnetic field are much smaller than the exchange coupling. This is manifest in (29) by the presence of the factor $1 - \kappa_\beta^2/14$ which implies that the contribution of spin misalignment may diverge when $\kappa_\beta$ is too large (i.e. for a strong field and/or large core anisotropy).

### 3. Magnetic SANS cross section

The quantity of interest in experimental SANS studies is the elastic magnetic differential scattering cross section $d\Sigma_M/d\Omega$, which is usually recorded on a two-dimensional position-sensitive detector. For the most commonly used scattering geometry in magnetic SANS experiments, where the applied magnetic field $\mathbf{B}_0 \parallel \mathbf{e}_z$ is perpendicular to the wavevector $\mathbf{k}_0 \parallel \mathbf{e}_x$ of the incident neutrons (see Fig. 4), $d\Sigma_M/d\Omega$ (for unpolarized neutrons) can be written as (Mühlbauer et al., 2019)

$$ \frac{d\Sigma_M}{d\Omega}(\mathbf{q}) = \frac{8\pi^3}{V} b_H^2 \left[ |\tilde{M}_x|^2 + |\tilde{M}_y|^2 \cos^2\theta_q + |\tilde{M}_z|^2 \sin^2\theta_q - (\tilde{M}_y\tilde{M}_z^* + \tilde{M}_y^*\tilde{M}_z) \sin\theta_q \cos\theta_q \right], \quad (31) $$

where $V$ is the scattering volume and $b_H = 2.91 \times 10^8$ Å$^{-1}$ m$^{-1}$ is the magnetic scattering length in the small-angle regime (the atomic magnetic form factor is approximated by 1,

<figure><img src="image_10.png" /><figcaption>Figure 4 A sketch of the perpendicular scattering geometry ($\mathbf{B}_0 \perp \mathbf{k}_0$). The scattering vector $\mathbf{q}$ corresponds to the difference between the wavevectors of the incident ($\mathbf{k}_0$) and scattered ($\mathbf{k}_1$) neutrons. The angle $\theta_q$ specifies the orientation of $\mathbf{q}$ on the detector. In the small-angle approximation, the component of $\mathbf{q}$ along $\mathbf{k}_0$ is neglected.</figcaption></figure>
since we are dealing with forward scattering). $\tilde{\mathbf{M}}(\mathbf{q}) = [\tilde{\mathbf{M}}_x(\mathbf{q}), \tilde{\mathbf{M}}_y(\mathbf{q}), \tilde{\mathbf{M}}_z(\mathbf{q})]$ represents the magnetization vector field $\mathbf{M}(\mathbf{r})$ in Fourier space, $\theta_q$ denotes the angle between the scattering vector $\mathbf{q}$ and $\mathbf{B}_0$ (not to be confused with the polar angle $\theta$ defined above), and the asterisk * stands for the complex conjugate. Note that in the perpendicular scattering geometry the Fourier components are evaluated in the plane $q_x = 0$.

The Fourier transform of the three-dimensional magnetization vector field (with a tilde above the symbol) is defined as

$$ \tilde{\mathbf{M}}(\mathbf{q}) = \frac{1}{(2\pi)^{3/2}} \int_V \mathbf{M}(\mathbf{r}) \exp(-i\mathbf{q} \cdot \mathbf{r}) \mathrm{d}^3\mathbf{r}, \quad (32) $$

$$ \mathbf{M}(\mathbf{r}) = \frac{1}{(2\pi)^{3/2}} \int_V \tilde{\mathbf{M}}(\mathbf{q}) \exp(i\mathbf{q} \cdot \mathbf{r}) \mathrm{d}^3\mathbf{q}. \quad (33) $$

For subsequent calculations, we introduce the following dimensionless quantities:

$$ \nu = \mathbf{q}R, \quad \tilde{\mathbf{M}} = \frac{(2\pi)^{3/2}}{4\pi R^3 M_0} \tilde{\mathbf{M}}, \quad (34) $$

and we express the dimensionless scattering vector in spherical coordinates as

$$ \mathbf{v} = [\nu \sin \theta_q \cos \phi_q, \nu \sin \theta_q \sin \phi_q, \nu \cos \theta_q]. \quad (35) $$

Next, in (32) we use the following first-order approximation for the real-space magnetization vector $\mathbf{m}(\xi)$ [see (5) and (7)]:

$$ \mathbf{m}(\xi) = \mathbf{m}_0 + \sum_{\beta=1}^{2} \mathbf{g}_{\beta} \psi_{\beta}(\xi). \quad (36) $$

As shown in Appendix C, the final expression for the Fourier transform of the magnetization is then given by

$$ \tilde{\mathbf{M}}(\nu) = \frac{j_1(\nu)}{\nu} \mathbf{m}_0 + \sum_{\beta=1}^{2} \mathbf{g}_{\beta} \sum_{\nu=1}^{\infty} \sum_{\mu=0}^{\nu} (-1)^{\nu} a_{\nu\mu}^{\beta} \rho_{\nu}^{\beta}(\nu) P_{2\nu}^{2\mu} (\cos \theta_q) \cos(2\mu\phi_q), \quad (37) $$

where

$$ \rho_{\nu}^{\beta}(\nu) = - \frac{\nu j_{2\nu-1}(\nu) \Upsilon_{\nu}(\kappa_{\beta}) - \kappa_{\beta} \Upsilon_{\nu}(\kappa_{\beta}) j_{2\nu}(\nu)}{\nu^2 + \kappa_{\beta}^2}, \quad (38) $$

$$ \mathfrak{N}_{\nu}(\tau) = \frac{\pi^{1/2}}{2} \sum_{s=0}^{\infty} \frac{(-1)^{\nu}(\tau/2)^{2(s+\nu)-1}}{s!\Gamma(2\nu + s + 1/2)} \quad (39) $$

and $\Upsilon_{\nu}(\kappa_{\beta})$ is given by (17). The zero-order term $\propto j_1(\nu)/\nu$ in (37) represents the form factor of a homogeneously magnetized sphere (Michels, 2021). In the limiting case of an infinite applied magnetic field, which is equivalent to the limit $\kappa_{\beta} \to \infty$, the additional terms [second line in (37)] vanish [compare with (25)] and the spherical form factor remains. On the other hand, if $k_s = 0$, the additional terms also vanish because, from the physical point of view, the Néel surface anisotropy cancels and from (18) we know that the coefficients $a_{\nu\mu}^{\beta}$ are linear in

$k_s$. Taking only the terms with $\nu = 1$ into account and setting $\phi_q = \pi/2$ ($\nu_x = 0$), corresponding to the scattering geometry where the applied magnetic field $\mathbf{B}_0 \parallel \mathbf{e}_z$ is perpendicular to the wavevector $\mathbf{k}_0 \parallel \mathbf{e}_x$ of the incident neutrons (Fig. 4), the expression for $\tilde{\mathbf{M}}(\nu)$ can be written as [compare with (21)]

$$ \tilde{\mathbf{M}}(\nu) \approx \frac{j_1(\nu)}{\nu} \mathbf{m}_0 - \frac{15k_s}{32} \sum_{\beta=1}^{2} \mathcal{R}_{\beta}(\nu) [\mathbf{g}_{\beta} \cdot \mathbf{V}(\theta_q, \pi/2) \cdot \mathbf{m}_0] \mathbf{g}_{\beta}, \quad (40) $$

where the radial function is

$$ \mathcal{R}_{\beta}(\nu) = \frac{\nu j_1(\nu)}{\nu^2 + \kappa_{\beta}^2} \frac{\Upsilon_1(\kappa_{\beta})}{\kappa_{\beta} \Upsilon'_1(\kappa_{\beta})} - \frac{j_2(\nu)}{\nu^2 + \kappa_{\beta}^2} \frac{\mathfrak{N}_1(\kappa_{\beta})}{\Upsilon'_1(\kappa_{\beta})}. \quad (41) $$

$\mathcal{R}_{\beta}(\nu)$ can be approximated for small $\kappa_{\beta}$ and, when only terms up to $s = 1$ in the infinite series (17) and (39) are kept,

$$ \mathcal{R}_{\beta}(\nu) = \frac{1}{4} \frac{\nu j_1(\nu)}{\nu^2 + \kappa_{\beta}^2} \frac{\kappa_{\beta}^2 + 14}{\kappa_{\beta}^2 + 7} - \frac{7}{4} \frac{j_2(\nu)}{\nu^2 + \kappa_{\beta}^2} \frac{\kappa_{\beta}^2 + 10}{\kappa_{\beta}^2 + 7}. \quad (42) $$

For small $\nu$ values, one finds the limit

$$ \lim_{\nu \to 0} \mathcal{R}_{\beta}(\nu) = 0, \quad (43) $$

which is consistent with

$$ \int_V \boldsymbol{\psi}(\mathbf{r}) \mathrm{d}^3\mathbf{r} = 0. \quad (44) $$

This can be seen by inspecting the definition of the Fourier transform in (32). Note that for $q \to 0$ the Fourier transform is proportional to the average of the magnetization vector field $\mathbf{M}$ and the maximum of this average is given by the homogeneous magnetization state. Using this result, the $\nu \to 0$ limit for the first-order approximation in $\boldsymbol{\psi}$ of the Fourier transform of the magnetization yields

$$ \lim_{\nu \to 0} \tilde{\mathbf{M}}(\nu, \theta_q, \phi_q) = \frac{1}{3} \mathbf{m}_0. \quad (45) $$

Beyond the linear approximation in $\boldsymbol{\psi}$, a non-vanishing term appears in $\tilde{\mathbf{M}}$ in the limit $\nu \to 0$, which reduces the Fourier components relative to the homogeneous magnetization state. In the second order in $\boldsymbol{\psi}$, the result is [compare (5)]

$$ \lim_{\nu \to 0} \tilde{\mathbf{M}}(\nu, \theta_q, \phi_q) = \left[ \frac{1}{3} - \frac{1}{2} \int_V \|\boldsymbol{\psi}(\xi)\|^2 \mathrm{d}^3\xi \right] \mathbf{m}_0. \quad (46) $$

Using (34) and

$$ \frac{\mathrm{d}\Sigma_M}{\mathrm{d}\Omega} = \frac{16\pi^2 R^6 M_0^2 b_H^2}{V} S_M, \quad (47) $$

the dimensionless two-dimensional magnetic SANS cross section $S_M(\nu, \theta_q)$ can be straightforwardly obtained as [compare (31)]

$$ S_M(\nu, \theta_q) = |\tilde{\mathbf{M}}_x|^2 + |\tilde{\mathbf{M}}_y|^2 \cos^2 \theta_q + |\tilde{\mathbf{M}}_z|^2 \sin^2 \theta_q - (\tilde{\mathbf{M}}_y \tilde{\mathbf{M}}_z^* + \tilde{\mathbf{M}}_y^* \tilde{\mathbf{M}}_z) \sin \theta_q \cos \theta_q. \quad (48) $$

In the limit $k_s \to 0$, the resulting cross section from (37) is
$$ \lim_{k_s \to 0} S_M(\nu, \theta_q) = \left[ \frac{j_1(\nu)}{\nu} \right]^2 (m_{0,x}^2 + m_{0,y}^2 \cos^2 \theta_q + m_{0,z}^2 \sin^2 \theta_q - 2m_{0,y}m_{0,z} \sin \theta_q \cos \theta_q). \quad (49) $$

Relation (49) nicely demonstrates that, depending on the orientation of the uniformly magnetized particle, different angular anisotropies become visible on the detector. For $\mathbf{m}_0 \parallel \mathbf{e}_x$ (i.e. $m_{0y} = m_{0z} = 0$) the scattering pattern is isotropic, while it exhibits a $\cos^2 \theta_q$ ($\sin^2 \theta_q$) type shape when $\mathbf{m}_0 \parallel \mathbf{e}_y$ ($\mathbf{m}_0 \parallel \mathbf{e}_z$).

Fig. 5 shows $S_M(\nu, \theta_q)$ along with the contribution of the individual Fourier components to (48). The upper row in Fig. 5 presents the results taking into account only the zero-order term $[j_1(\nu)/\nu]\mathbf{m}_0$ from (40), while in the lower row the second-order term ($\nu = 1$) is additionally included. Since the zero-order term represents the case of a homogeneously magnetized NM, this comparison provides useful insights about the impact of the Néel surface anisotropy on the magnetic SANS cross section. In the case of a uniformly magnetized NM (upper row) the Fourier components $|\tilde{M}_x|^2$, $|\tilde{M}_y|^2$ and $|\tilde{M}_z|^2$ are isotropic (rotational symmetry), while including the second-order terms (lower row) leads to anisotropic behaviour of the transverse components $|\tilde{M}_x|^2$ and $|\tilde{M}_y|^2$. The cross term (CT) averages to zero for both situations and the dominant contribution to the magnetic SANS cross section (for the parameters chosen in Fig. 5) is given by the $|\tilde{M}_z|^2$ component. Therefore, it may be concluded that the impact of the Néel surface anisotropy on $S_M(\nu, \theta_q)$ is relatively small. By comparing the $S_M(\nu, \theta_q)$ from the upper and lower rows, it is seen that by including the Néel surface anisotropy the circular symmetry of the zeros of $S_M$ (deep-blue colours) is broken. This feature becomes more clearly visible by analyzing the azimuthal average of $S_M(\nu, \theta_q)$, which is readily computed as

$$ \mathcal{I}(\nu) = \frac{1}{2\pi} \int_{0}^{2\pi} S_M(\nu, \theta_q) \, d\theta_q. \quad (50) $$

In the limit $k_s \to 0$, the azimuthal average corresponding to (49) is

$$ \lim_{k_s \to 0} \mathcal{I}(\nu) = \left[ \frac{j_1(\nu)}{\nu} \right]^2 \frac{\|\mathbf{m}_0\|^2 + (\mathbf{e}_x \cdot \mathbf{m}_0)^2}{2}. \quad (51) $$

We have also calculated the pair-distance distribution function

$$ \mathcal{P}(\xi) = \xi^2 \int_0^\infty \mathcal{I}(\nu) j_0(\nu\xi) \, \nu^2 \, d\nu \quad (52) $$

and the correlation function

$$ \mathcal{C}(\xi) = \mathcal{P}(\xi)/\xi^2. \quad (53) $$

In the limit $k_s \to 0$, the pair-distance distribution and the correlation function corresponding to (51) are

$$ \lim_{k_s \to 0} \mathcal{P}(\xi) = \frac{\pi \xi^2}{6} \left(1 - \frac{3\xi}{4} + \frac{\xi^3}{16}\right) \frac{\|\mathbf{m}_0\|^2 + (\mathbf{e}_x \cdot \mathbf{m}_0)^2}{2}, \quad (54) $$

$$ \lim_{k_s \to 0} \mathcal{C}(\xi) = \frac{\pi}{6} \left(1 - \frac{3\xi}{4} + \frac{\xi^3}{16}\right) \frac{\|\mathbf{m}_0\|^2 + (\mathbf{e}_x \cdot \mathbf{m}_0)^2}{2}. \quad (55) $$

These functions are displayed graphically in Fig. 6. Due to the surface-anisotropy-induced spin disorder, the form-factor extrema of $\mathcal{I}(\nu)$ [Fig. 6(a)] are damped and shifted slightly to larger $q$ values [i.e. smaller structures; compare the first minimum of $\mathcal{I}(\nu)$ for $k_s = 3$]. Moreover, as already observed in numerical micromagnetic continuum simulations (Vivas et al., 2017, 2020), the oscillations are damped for the case of surface spin disorder, which mimics the effect of a particle-size distribution or of instrumental resolution. In agreement with

<figure><img src="image_11.png" /><figcaption>Figure 5</figcaption></figure>

Results for the two-dimensional Fourier components $|\tilde{M}_x|^2$, $|\tilde{M}_y|^2$, $|\tilde{M}_z|^2$, CT = $-(\tilde{M}_y, \tilde{M}_z^* + \tilde{M}_y^*, \tilde{M}_z)$ and for the total magnetic SANS cross section $S_M(\nu, \theta_q)$ [(48)] using expression (40). The upper row shows the results taking into account only the zero-order term in (40), which corresponds to the case of a homogeneously magnetized particle. The lower row displays the results when the second-order term ($\nu = 1$) in (40) is taken into account. The parameters are $\mathbf{e}_x = \mathbf{e}_z$, $\mathbf{b}_0 = 0.1\mathbf{e}_z$ ($B_0 \simeq 48$ mT), $k_c = 0.1$, $k_s = 3$ and $\mathbf{m}_0 = [\sin \beta \cos \alpha, \sin \beta \sin \alpha, \cos \beta]$. Note that $\nu_y$ and $\nu_z$ denote the dimensionless components of the scattering vector [compare equation (34)]. Since the Néel surface anisotropy effectively has a cubic symmetry (see Fig. 1), we average $S_M$ over the angles $\alpha = (45^\circ, 135^\circ, 225^\circ, 315^\circ)$ and $\beta = 20^\circ$. A logarithmic colour scale is used.
this observation is the finding that the maximum of the $P(\xi)$ function [Fig. 6(b)] appears at smaller distances $\xi$ than in the homogeneous case. Likewise, due to spin disorder, the $C(\xi)$ function [Fig. 6(c)] exhibits a larger amplitude (Mettus & Michels, 2015).

To analyze the role of the surface anisotropy more quantitatively, we have computed the following quantities, which describe the deviation of the one-dimensional SANS cross section and of the pair-distance distribution function from the homogeneous particle case:

$$ \epsilon\mathcal{I}(k_s) = \frac{\int_0^\infty |\mathcal{I}(k_s = 0, \nu) - \mathcal{I}(k_s, \nu)| d\nu}{\int_0^\infty |\mathcal{I}(k_s = 0, \nu)| d\nu}, \quad (56) $$

$$ \epsilon\mathcal{P}(k_s) = \frac{\int_0^2 |\mathcal{P}(k_s = 0, \xi) - \mathcal{P}(k_s, \xi)| d\xi}{\int_0^2 |\mathcal{P}(k_s = 0, \xi)| d\xi}. \quad (57) $$

Fig. 6(d) depicts both $\epsilon\mathcal{I}(k_s)$ and $\epsilon\mathcal{P}(k_s)$ as a function of $k_s$. The difference is only of the order of a few percent, which suggests that the effect of surface anisotropy on the SANS observables is relatively weak within the present analytical approximation: see our accompanying numerical work (Adams et al., 2022), which takes into account the full nonlinearity of the micromagnetic equations. However, this is only true for the magnetic interactions considered here. Taking into account the anisotropic and long-range dipole-dipole interaction and the asymmetric Dzyaloshinskii-Moriya interaction will very likely result in more inhomogeneous spin structures

<figure><img src="image_12.png" /><figcaption>Figure 6 Results for (a) the azimuthally averaged SANS cross section $\mathcal{I}(\nu)$, (b) the pair-distance distribution function $\mathcal{P}(\xi)$, (c) the correlation function $\mathcal{C}(\xi)$, and (d) the quantities $\epsilon\mathcal{I}(k_s)$ [(56)] and $\epsilon\mathcal{P}(k_s)$ [(57)] ($B_0 \approx 48$ mT). For the homogeneous case [blue curves in panels (a)-(c)] the surface anisotropy is set to $k_s = 0$, and for the inhomogeneous case [red curves in panels (a)-(c)] we use the same parameters as in Fig. 5. The functional dependence of $\mathcal{I}(\nu)$, $\mathcal{P}(\xi)$ and $\mathcal{C}(\xi)$ for the uniformly magnetized particle all correspond to the analytically well known cases, i.e. (51), (54) and (55).</figcaption></figure>

and in larger deviations from the macrospin model (Vivas et al., 2017, 2020; Pathak & Hertel, 2021). Likewise, for NMs of elongated shapes, the surface anisotropy renders an additional first-order contribution to the effective energy (Garanin & Kachkachi, 2003), in addition to the second-order cubic contribution discussed above. This new shape-induced contribution could also lead to an enhancement of the spin misalignment. The analytical calculations presented here provide a general framework for future studies of more complicated (anisotropic) magnetic interactions. The approach can be straightforwardly adapted to other particle shapes such as a circular planar disc.

## 4. Conclusions

We have analytically computed the magnetization distribution and the ensuing magnetic SANS cross section of a spherical nanoparticle. Our micromagnetic Hamiltonian takes into account the isotropic exchange interaction, an external magnetic field, a uniaxial anisotropy for the particle's core and Néel anisotropy on its boundary. The resulting Helmholtz equation has been solved by expanding the real-space magnetization in terms of spherical Bessel functions and spherical harmonics. The central results are the infinite series (16) and its second-order expansion (21) for the real-space magnetization, and the corresponding Fourier transforms (37) and (40). Using these expressions, the two-dimensional magnetic SANS cross section $S_M(\nu, \theta_q)$, the azimuthally averaged SANS signal $\mathcal{I}(\nu)$, and the correlation functions $\mathcal{P}(\xi)$ and $\mathcal{C}(\xi)$ have been obtained and compared with the case of a homogeneous spin configuration (uniform magnetization vector field). The signature of Néel surface anisotropy (of constant $k_s$) has been identified in all of these functions. However, its effect is relatively small, even for large values of $k_s$. Taking into account the magnetodipolar and/or the Dzyaloshinskii-Moriya interaction, or shape asymmetry, will probably result in configurations with stronger spin misalignment (e.g. in vortex-type textures or skyrmions) and thereby in more prominent signatures in the SANS cross section and correlation function. These interactions are beyond the scope of the current analytical approach and will be considered in our future (numerical) work.

# APPENDIX A

## Solution of the boundary value problem of the Helmholtz equation

The coefficients $c_{\ell m}^{\beta}$ in the fundamental solution (15) of the Helmholtz equation (9) must be determined such that the Neumann boundary condition (10) is satisfied. For this purpose, we use the method of least squares, where we make use of the orthogonality properties of the spherical harmonics $Y_{\ell m}(\theta, \phi)$. The normal derivative of (15) at the surface of the NM ($\xi = 1$) is
$$ \left. \frac{d\psi_{\beta}}{d\xi} \right|_{\xi=1} = \sum_{\ell=0}^{\infty} \sum_{m=-\ell}^{\ell} c_{\ell m}^{\beta} [j_{\ell}(i\kappa_{\beta}\xi)]'_{\xi=1} Y_{\ell m}(\theta, \phi), \quad (58) $$

where

$$ [j_{\ell}(i\kappa_{\beta}\xi)]'_{\xi=1} = \left. \frac{d}{d\xi} j_{\ell}(i\kappa_{\beta}\xi) \right|_{\xi=1}. \quad (59) $$

Our goal is now to minimize the following error functional with respect to the coefficients $c_{\ell m}^{\beta}$:

$$ \begin{aligned} \epsilon[c_{\ell m}^{\beta}] = & \int_{0}^{2\pi} \int_{0}^{\pi} \left| \sum_{\alpha \in \{x, y, z\}} \chi_{\alpha}^{\beta} |n_{\alpha}| \right. \\ & \left. - \sum_{\ell=0}^{\infty} \sum_{m=-\ell}^{\ell} c_{\ell m}^{\beta} [j_{\ell}(i\kappa_{\beta}\xi)]'_{\xi=1} Y_{\ell m}(\theta, \phi) \right|^2 \sin \theta \, d\theta \, d\phi. \end{aligned} \quad (60) $$

The minimum of this error functional is found from the condition that the partial derivatives of $\epsilon[c_{\ell m}^{\beta}]$ with respect to $(c_{ij}^{\beta})^*$ vanish,

$$ \frac{\partial \epsilon[c_{\ell m}^{\beta}]}{\partial (c_{ij}^{\beta})^*} = 0, \quad (61) $$

where * denotes the complex conjugate. Using the orthogonality relation (62) [p. 378 (14.30.8) of Olver et al. (2010)] and by defining the integral $I_{ij}^{\alpha}$ by (63),

$$ \int_{0}^{2\pi} \int_{0}^{\pi} Y_{\ell m}(\theta, \phi) Y_{ij}^{*}(\theta, \phi) \sin \theta \, d\theta \, d\phi = \delta_{\ell l} \delta_{mj}, \quad (62) $$

$$ \int_{0}^{2\pi} \int_{0}^{\pi} Y_{ij}^{*}(\theta, \phi) |n_{\alpha}| \sin \theta \, d\theta \, d\phi = I_{ij}^{\alpha}, \quad (63) $$

the solution of (61) can be written as

$$ c_{ij}^{\beta} = \frac{1}{[j_{\ell}(i\kappa_{\beta}\xi)]'_{\xi=1}} \sum_{\alpha \in \{x, y, z\}} \chi_{\alpha}^{\beta} I_{ij}^{\alpha}. \quad (64) $$

Alternatively, again using the indices $\ell$ and $m$, expressing the coefficient $\chi_{\alpha}^{\beta}$ as in (13) and using the matrix vector product, the coefficients $c_{\ell m}^{\beta}$ can be more explicitly written as

$$ c_{\ell m}^{\beta} = \frac{k_s \mathbf{g}_{\beta} \cdot \mathrm{diag}[I_{\ell m}^{x}, I_{\ell m}^{y}, I_{\ell m}^{z}] \cdot \mathbf{m}_0}{[j_{\ell}(i\kappa_{\beta}\xi)]'_{\xi=1}}. \quad (65) $$

For some low orders of $\ell$ and $m$, the exact solutions of the integrals $I_{\ell m}^{\alpha}$ are presented in Appendix B.

From (65) several conclusions can be drawn. First, the zero-order term ($\ell = m = 0$) in (15) vanishes, which can easily be shown by rewriting the coefficient $c_{00}^{\beta}$ as

$$ c_{00}^{\beta} = \frac{k_s(\pi^{1/2})(\mathbf{g}_{\beta} \cdot \mathbf{m}_0)}{[j_0(i\kappa_{\beta}\xi)]'_{\xi=1}} = 0, \quad (66) $$

and this is due to the orthogonality of $\mathbf{m}_0$ and $\mathbf{g}_{\beta}$ with $\beta \in \{1, 2\}$ [see (6)]. Moreover, we see that

$$ \psi_{\beta}(\xi = 0, \theta, \phi) = 0, \quad (67) $$

which is a consequence of the behaviour of the spherical Bessel functions of the first kind at the origin. Since the $j_{\ell}$ with $\ell \ge 1$ all vanish at the origin $\xi = 0$, this implies that the total $\psi_{\beta}$ also vanishes at the origin. Note that $j_0(0) = 1$ but does not

contribute to $\psi_{\beta}$ since $c_{00}^{\beta} = 0$. From the physical point of view this makes sense, since the spin misalignment is caused by the Néel surface anisotropy and thus, for symmetry reasons, there is no spin disorder at the centre of the spherical NM; the highest misalignment is found at its surface. Secondly, in Table 1 in Appendix B it is seen that the coefficients $I_{\ell m}^{\alpha}$ vanish for odd $\ell$ or $m$ and that $I_{\ell m}^{\alpha} = I_{\ell, -m}^{\alpha}$, so that the expansion coefficients also exhibit this symmetry,

$$ c_{\ell m}^{\beta} = c_{\ell, -m}^{\beta}. \quad (68) $$

Taking these properties into account, one can express the solution (15) more conveniently in terms of the associated Legendre polynomials $P_{\ell m}^{m}(\cos \theta)$ with $\ell = 2\nu$ and $m = 2\mu$ [note that we use the convention that $Y_{\ell m}(\theta, \phi) = N_{\ell m} P_{\ell}^{m}(\cos \theta) \times \exp(i\mu\phi)$ [p. 378 (14.30.1) of Olver et al. (2010)]],

$$ \psi_{\beta} = \sum_{\nu=1}^{\infty} \sum_{\mu=0}^{\nu} a_{\nu\mu}^{\beta} \Upsilon_{\nu}(\kappa_{\beta}\xi) P_{2\nu}^{2\mu}(\cos \theta) \cos(2\mu\phi), \quad (69) $$

where $\Upsilon_{\nu}(\tau)$ is defined in (17) and $a_{\nu\mu}^{\beta}$ in (18).

The infinite series in (18) is a consequence of the relation between the spherical Bessel functions of the first kind and the ordinary Bessel functions of the first kind [p. 262 (10.47.3) of Olver et al. (2010)],

$$ j_{\ell}(\tau) = \left(\frac{\pi}{2\tau}\right)^{1/2} J_{l+1/2}(\tau), \quad (70) $$

and the well known series representation [p. 262 (10.47.3) of Olver et al. (2010)]

$$ J_{l}(\tau) = \sum_{s=0}^{\infty} \frac{(-1)^{s}(\tau/2)^{2s+l}}{s!\Gamma(l+s+1)}. \quad (71) $$

## APPENDIX B: Integral coefficients

The integrals $I_{\ell m}^{\alpha}$ can be simplified, since both $n_{\alpha}$ and $Y_{\ell m}^{*}(\theta, \phi)$ are separable functions in $\theta$ and $\phi$. By rewriting the spherical harmonics in terms of complex exponentials and associated Legendre polynomials, and by using the definition of the Cartesian components of the surface normal vector $\mathbf{n}$ from (3) [note that we use the convention that $Y_{\ell m}(\theta, \phi) = N_{\ell m} P_{\ell}^{m}(\cos \theta) \exp(i\mu\phi)$ [p. 378 (14.30.1) of Olver et al. (2010)]], the integrals are expressed as follows:

$$ I_{\ell m}^{x} = N_{\ell m} K_{\ell m}^{(1)} U_{m}^{(1)}, \quad (72) $$

$$ I_{\ell m}^{y} = N_{\ell m} K_{\ell m}^{(1)} U_{m}^{(2)}, \quad (73) $$

$$ I_{\ell m}^{z} = N_{\ell m} K_{\ell m}^{(2)} U_{m}^{(3)}, \quad (74) $$

where

$$ K_{\ell m}^{(1)} = \int_{0}^{\pi} P_{\ell}^{m}(\cos \theta) |\sin \theta| \sin \theta \, d\theta, \quad (75) $$
$$K_{\ell m}^{(2)} = \int_{0}^{\pi} P_{\ell}^{m}(\cos \theta) |\cos \theta| \sin \theta \, d\theta, \quad (76)$$

$$U_m^{(1)} = \int_0^{2\pi} \exp(-im\phi) |\cos \phi| \, d\phi = \frac{2[1 + (-1)^{|m|}]}{1 - m^2 + \delta_{|m|,1}}, \quad (77)$$

$$U_m^{(2)} = \int_0^{2\pi} \exp(-im\phi) |\sin \phi| \, d\phi = \frac{2[1 + (-1)^{|m|}]}{1 - m^2 + \delta_{|m|,1}}, \quad (78)$$

$$U_m^{(3)} = \int_0^{2\pi} \exp(-im\phi) \, d\phi = 2\pi\delta_{0,m}, \quad (79)$$

$$N_{\ell m} = \left[ \frac{2\ell + 1}{4\pi} \frac{(\ell - m)!}{(\ell + m)!} \right]^{1/2}. \quad (80)$$

The integrals $U_m^{(1)}$, $U_m^{(2)}$ and $U_m^{(3)}$ are solvable straightforwardly using Euler's formula for the complex exponential and by splitting the region of integration according to the absolute values of the trigonometric functions. In the denominator of $U_m^{(1)}$ and $U_m^{(2)}$ we have included the Kronecker delta $\delta_{|m|,1}$ to take account of the cases $m = \pm 1$. It is common to express the integrals $K_{\ell m}^{(1)}$ and $K_{\ell m}^{(2)}$ by the substitution $x = \cos \theta$ ($dx = -\sin \theta \, d\theta$),

$$K_{\ell m}^{(1)} = \int_{-1}^{1} P_{\ell}^{m}(x) (1 - x^2)^{1/2} \, dx, \quad (81)$$

$$K_{\ell m}^{(2)} = \int_{-1}^{1} P_{\ell}^{m}(x) |x| \, dx. \quad (82)$$

Using $U_m^{(3)} = 2\pi\delta_{0,m}$, we need only compute $K_{\ell,0}^{(2)}$ such that the associated Legendre polynomials in $K_{\ell m}^{(2)}$ are reduced to the Legendre polynomials with one index only [p. 352 of Olver et al. (2010)]. By considering the symmetry properties of the Legendre polynomials it becomes clear that the integral must vanish for odd $\ell$ and can be simplified for even $\ell$ in the following way:

$$K_{\ell,0}^{(2)} = [1 + (-1)^{\ell}] \int_{0}^{1} P_{\ell}(x) x \, dx. \quad (83)$$

The closed form of $K_{\ell,0}^{(2)}$ is then found in terms of the Gamma function [p. 771 (7.126.1) of Gradshteyn & Ryzhik (2007)],

$$K_{\ell,0}^{(2)} = \frac{(\pi^{1/2})[1 + (-1)^{\ell}]}{4\Gamma(3/2 - \ell/2)\Gamma(\ell/2 + 2)}. \quad (84)$$

The overall solution for $I_{\ell m}^z$ is then written as

$$I_{\ell m}^z = \frac{\pi[1 + (-1)^{\ell}] (2\ell + 1)^{1/2}}{4\Gamma(3/2 - \ell/2)\Gamma(\ell/2 + 2)} \delta_{0,m}. \quad (85)$$

Since we only have to calculate $K_{\ell m}^{(1)}$ for even $m$ [because $U_m^{(1)}$ and $U_m^{(2)}$ include the term $1 + (-1)^{|m|}$], we may use the index $m = 2\mu$ and, by exploiting the parity of the associated Legendre polynomials $P_{\ell}^{m}(-x) = (-1)^{\ell+m}P_{\ell}^{m}(x)$, we find

<table><caption>**Table 1**

Values of the integrals (19) for some small values of $\ell$ and $m$.</caption><tbody><tr><td colspan="8">I<sup>*</sup><sub>ℓm</sub> = ∫<sub>0</sub><sup>2π</sup> ∫<sub>0</sub><sup>π</sup> Y<sup>*</sup><sub>ℓm</sub>(θ, φ) |sin θ cos φ| sin θ dθ dφ = (-1)<sup>|m|/2</sup> × <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">2[1 + (-1)<sup>|m|</sup>][1 + (-1)<sup>ℓ</sup>]</div><div style="padding: 0 5px;">4π</div></div> <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">2ℓ + 1</div><div style="padding: 0 5px;">(ℓ - m)!</div></div> <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">1 - m<sup>2</sup> + δ<sub>|m|,1</sub></div><div style="padding: 0 5px;">4π</div></div> ∫<sub>0</sub><sup>1</sup> P<sup>m</sup><sub>ℓ</sub>(x) (1 - x<sup>2</sup>)<sup>1/2</sup> dx</td></tr><tr><td colspan="8">m</td></tr><tr><td>ℓ</td><td>0</td><td>±1</td><td>±2</td><td>±3</td><td>±4</td><td>±5</td><td>±6</td></tr><tr><td>0</td><td colspan="7">π<sup>1/2</sup></td></tr><tr><td>1</td><td colspan="7">0</td></tr><tr><td>2</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">(5π)<sup>1/2</sup></div><div style="padding: 0 5px;">8</div></div></td><td>0</td><td>0</td><td>1/8 (15π/2)<sup>1/2</sup></td><td colspan="4"></td></tr><tr><td>3</td><td colspan="7">0</td></tr><tr><td>4</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">3(π<sup>1/2</sup>)</div><div style="padding: 0 5px;">64</div></div></td><td>0</td><td>0</td><td>1/32 (5π/2)<sup>1/2</sup></td><td>0</td><td>-1/64 (35π/2)<sup>1/2</sup></td><td colspan="2"></td></tr><tr><td>5</td><td colspan="7">0</td></tr><tr><td>6</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">5(13π)<sup>1/2</sup></div><div style="padding: 0 5px;">1024</div></div></td><td>0</td><td>0</td><td>(1365π)<sup>1/2</sup>/2048</td><td>0</td><td>-3/1024 (91π/2)<sup>1/2</sup></td><td>0</td><td>(3003π)<sup>1/2</sup>/2048</td></tr><tr><td colspan="8">I<sup>*</sup><sub>ℓm</sub> = ∫<sub>0</sub><sup>2π</sup> ∫<sub>0</sub><sup>π</sup> Y<sup>*</sup><sub>ℓm</sub>(θ, φ) |sin θ sin φ| sin θ dθ dφ = <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">2[1 + (-1)<sup>|m|</sup>][1 + (-1)<sup>ℓ</sup>]</div><div style="padding: 0 5px;">1 - m<sup>2</sup> + δ<sub>|m|,1</sub></div></div> <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">2ℓ + 1</div><div style="padding: 0 5px;">(ℓ - m)!</div></div> <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">1 - m<sup>2</sup> + δ<sub>|m|,1</sub></div><div style="padding: 0 5px;">4π</div></div> ∫<sub>0</sub><sup>1</sup> P<sup>m</sup><sub>ℓ</sub>(x) (1 - x<sup>2</sup>)<sup>1/2</sup> dx</td></tr><tr><td colspan="8">m</td></tr><tr><td>ℓ</td><td>0</td><td>±1</td><td>±2</td><td>±3</td><td>±4</td><td>±5</td><td>±6</td></tr><tr><td>0</td><td colspan="7">π<sup>1/2</sup></td></tr><tr><td>1</td><td colspan="7">0</td></tr><tr><td>2</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">(5π)<sup>1/2</sup></div><div style="padding: 0 5px;">8</div></div></td><td>0</td><td>0</td><td>1/8 (15π/2)<sup>1/2</sup></td><td colspan="4"></td></tr><tr><td>3</td><td colspan="7">0</td></tr><tr><td>4</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">3(π<sup>1/2</sup>)</div><div style="padding: 0 5px;">64</div></div></td><td>0</td><td>0</td><td>1/32 (5π/2)<sup>1/2</sup></td><td>0</td><td>-1/64 (35π/2)<sup>1/2</sup></td><td colspan="2"></td></tr><tr><td>5</td><td colspan="7">0</td></tr><tr><td>6</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">5(13π)<sup>1/2</sup></div><div style="padding: 0 5px;">1024</div></div></td><td>0</td><td>0</td><td>(1365π)<sup>1/2</sup>/2048</td><td>0</td><td>-3/1024 (91π/2)<sup>1/2</sup></td><td>0</td><td>(3003π)<sup>1/2</sup>/2048</td></tr><tr><td colspan="8">I<sup>*</sup><sub>ℓm</sub> = ∫<sub>0</sub><sup>2π</sup> ∫<sub>0</sub><sup>π</sup> Y<sup>*</sup><sub>ℓm</sub>(θ, φ) |cos θ| sin θ dθ dφ = <div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">π[1 + (-1)<sup>ℓ</sup>](2ℓ + 1)<sup>1/2</sup></div><div style="padding: 0 5px;">4Γ(3/2 - ℓ/2)Γ(ℓ/2 + 2)</div></div> δ<sub>m,0</sub></td></tr><tr><td colspan="8">m</td></tr><tr><td>ℓ</td><td>0</td><td>±1</td><td>±2</td><td>±3</td><td>±4</td><td>±5</td><td>±6</td></tr><tr><td>0</td><td colspan="7">π<sup>1/2</sup></td></tr><tr><td>1</td><td colspan="7">0</td></tr><tr><td>2</td><td><div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">(5π)<sup>1/2</sup></div><div style="padding: 0 5px;">4</div></div></td><td>0</td><td>0</td><td colspan="4"></td></tr><tr><td>3</td><td colspan="7">0</td></tr><tr><td>4</td><td>-<div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">π<sup>1/2</sup></div><div style="padding: 0 5px;">8</div></div></td><td>0</td><td>0</td><td>0</td><td>0</td><td colspan="2"></td></tr><tr><td>5</td><td colspan="7">0</td></tr><tr><td>6</td><td><div style="display: inline-block; vertical-align: middle; text-align: center;"><div style="border-bottom: 1px solid black; padding: 0 5px;">(13π)<sup>1/2</sup></div><div style="padding: 0 5px;">64</div></div></td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td><td>0</td></tr></tbody></table>
$$K_{\ell, 2\mu}^{(1)} = [1 + (-1)^{\ell}] \int_{0}^{1} P_{\ell}^{2\mu}(x) (1 - x^2)^{1/2} dx. \quad (86)$$

From these results it is seen that the integrals $I_{\ell m}^{\alpha}$ with $\alpha \in \{x, y, z\}$ vanish for odd $\ell$ and $m$ [note that in (86) this is only the case for $m = 2\mu$]. For the remaining integrals $K_{2\nu, 2\mu}^{(1)}$ in (86), where $\ell = 2\nu$, we do not give an expression in closed form. However, one must exist in terms of the Gamma function or the Beta function, since [p. 324 (3.251.2) of Gradshteyn & Ryzhik (2007)]

$$\int_{0}^{1} x^{\ell} (1 - x^{2})^{1/2} dx = \frac{1}{2} B\left(\frac{s+1}{2}, \frac{3}{2}\right), \quad (87)$$

where $B(\cdot, \cdot)$ is the Beta function (Euler integral), and the associated Legendre functions $P_{2\nu}^{2\mu}$ of even order and degree are true polynomials, as seen for example from the related Rodrigues formula [p. 360 (14.7.14) of Olver et al. (2010)]. We used Mathematica (Wolfram Inc.) to determine the integrals up to the sixth order in $\ell$ (see Table 1).

## APPENDIX C

### Derivation of the Fourier transform of the magnetization

The Fourier transform of the magnetization vector field $M(r)$ is written as

$$\tilde{M}(q) = \frac{1}{(2\pi)^{3/2}} \int_V M(r) \exp(-iq \cdot r) d^3 r. \quad (88)$$

In the following, we will use dimensionless quantities. For this purpose, we define the dimensionless scattering vector $v = qR$, where $R$ is the radius of the nanomagnet, the dimensionless position vector $\xi = r/R$ and the dimensionless magnetization vector $m = M/M_0$, where $M_0$ is the saturation magnetization. Substituting in (88) results in

$$\tilde{M}(v) = \frac{R^3 M_0}{(2\pi)^{3/2}} \int_V m(\xi) \exp(-iv \cdot \xi) d^3 \xi. \quad (89)$$

The dimensionless Fourier transform $\tilde{M}$ is then defined as

$$\tilde{M}(v) = \frac{1}{4\pi} \int_V m(\xi) \exp(-iv \cdot \xi) d^3 \xi, \quad (90)$$

with

$$\tilde{M} = \frac{4\pi R^3 M_0}{(2\pi)^{3/2}} \tilde{M}. \quad (91)$$

The next step consists of calculating the Fourier integral of the first-order approximation (36) of the magnetization vector $m$. Since (36) is formulated in dimensionless spherical coordinates $\xi, \theta, \phi$, it is convenient to express the scattering vector $v$ in spherical coordinates as well,

$$\xi = [\xi \sin \theta \cos \phi, \xi \sin \theta \sin \phi, \xi \cos \theta], \quad (92)$$

$$v = [\nu \sin \theta_q \cos \phi_q, \nu \sin \theta_q \sin \phi_q, \nu \cos \theta_q], \quad (93)$$

so that the plane-wave expansion of the complex exponential can be used (Jackson, 1999),

$$\exp(-iv \cdot \xi) = 4\pi \sum_{k=0}^{\infty} \sum_{n=-k}^{k} (-i)^k j_k(\nu\xi) Y_{kn}^*(\theta, \phi) Y_{kn}(\theta_q, \phi_q). \quad (94)$$

The Fourier integral (90) is then expanded into the following infinite series:

$$\tilde{M}(v) = \sum_{k=0}^{\infty} \sum_{n=-k}^{k} (-i)^k \hat{M}_{kn}(v) Y_{kn}(\theta_q, \phi_q), \quad (95)$$

where

$$\hat{M}_{kn}(v) = \int_0^{2\pi} \int_0^\pi \int_0^1 m(\xi) j_k(\nu\xi) Y_{kn}^*(\theta, \phi) \xi^2 \sin\theta \, d\xi \, d\theta \, d\phi. \quad (96)$$

We now use the infinite series (15) for $\psi_\beta$ to express the first-order approximation of the magnetization, which leads to

$$\mathbf{m}(\xi) = \mathbf{m}_0 + \sum_{\beta=1}^{2} \mathbf{g}_{\beta} \sum_{\ell=0}^{\infty} \sum_{m=-\ell}^{\ell} c_{\ell m}^{\beta} j_{\ell}(i\kappa_{\beta}\xi) Y_{\ell m}(\theta, \phi). \quad (97)$$

Since the integral transform (96) is linear, each term in (97) can be separately transformed. For the zero-order term, we obtain

$$\int_0^{2\pi} \int_0^\pi \int_0^1 \mathbf{m}_0 j_k(\nu\xi) Y_{kn}^*(\theta, \phi) \xi^2 \sin\theta \, d\xi \, d\theta \, d\phi \\ = (4\pi)^{1/2} \frac{j_1(\nu)}{\nu} \mathbf{m}_0 \delta_{k,0} \delta_{n,0}. \quad (98)$$

This result is well known to the neutron-scattering community as the spherical form factor, corresponding to a uniformly magnetized spherical particle (Michels, 2021). In the second step, we carry out the integration of the higher-order terms from (97). The radial and angular parts in the higher-order terms of (97) are multiplicative, such that the volume integral (96) is separable into

$$A_{\ell k}^{\beta}(v) = \int_{0}^{1} j_{\ell}(i\kappa_{\beta}\xi) j_{k}(\nu\xi) \xi^2 d\xi, \quad (99)$$

$$B_{\ell m}^{kn} = \int_{0}^{2\pi} \int_{0}^{\pi} Y_{\ell m}(\theta, \phi) Y_{kn}^{*}(\theta, \phi) \sin\theta \, d\theta \, d\phi. \quad (100)$$

As an intermediate result, the integral (96) is then rewritten as

$$\hat{M}_{kn}(v) = (4\pi)^{1/2} \frac{j_1(v)}{v} \mathbf{m}_0 \delta_{k,0} \delta_{n,0} \\ + \sum_{\beta=1}^{2} \mathbf{g}_{\beta} \sum_{\ell=0}^{\infty} \sum_{m=-\ell}^{\ell} c_{\ell m}^{\beta} A_{\ell k}^{\beta}(v) B_{\ell m}^{kn}. \quad (101)$$

The integral (100) directly corresponds to the orthogonality relation of the spherical harmonics and thereby we have $B_{\ell m}^{kn} = \delta_{\ell k} \delta_{mn}$ [p. 378 (14.30.8) of Olver et al. (2010)]. Due to the term $\delta_{\ell k}$ in $B_{\ell m}^{kn}$ and since $B_{\ell m}^{kn}$ and $A_{\ell k}^{\beta}$ are multiplicative in (101), the following spherical Hankel transform results:
$$A_{kk}^{\beta}(\nu) = \int_{0}^{1} j_k(i\kappa_{\beta}\xi) j_k(\nu\xi) \xi^2 d\xi. \qquad (102)$$

In order to calculate these integrals, we replace [using (70)] the spherical Bessel functions of the first kind $j_k(\cdot)$ with the ordinary Bessel functions of the first kind $J_l(\cdot)$. This yields

$$A_{kk}^{\beta}(\nu) = \left(\frac{\pi}{2\nu}\right)^{1/2} \left(\frac{\pi}{2i\kappa_{\beta}}\right)^{1/2} W_k^{\beta}(\nu), \qquad (103)$$

where

$$W_k^{\beta}(\nu) = \int_0^1 j_{k+1/2}(i\kappa_{\beta}\xi) J_{k+1/2}(\nu\xi) \xi \, d\xi \qquad (104)$$

is the Hellom integral. The result for the integral $W_k^\beta$ can be found in the textbook by Gradshteyn & Ryzhik (2007) [p. 664 (6.521.1)] and is written as

$$W_k^{\beta}(\nu) = - \frac{\nu J_{k-1/2}(\nu) J_{k+1/2}(i\kappa_{\beta}) - i\kappa_{\beta} J_{k-1/2}(i\kappa_{\beta}) J_{k+1/2}(\nu)}{\nu^2 + \kappa_{\beta}^2}. \qquad (105)$$

Again, upon applying relation (70), we write the solution for $A_{kk}^\beta(\nu)$ more conveniently in terms of spherical Bessel functions of the first kind as

$$A_{kk}^{\beta}(\nu) = - \frac{\nu j_{k-1}(\nu) j_k(i\kappa_{\beta}) - i\kappa_{\beta} j_{k-1}(i\kappa_{\beta}) j_k(\nu)}{\nu^2 + \kappa_{\beta}^2}. \qquad (106)$$

Now that the integrals have been obtained, we have to substitute (101) in (95). The term on the first line of (101) only accounts for $k=0$ and $n=0$ in (95). Therefore, the contribution of this term to $\tilde{M}$ is $[j_1(\nu)/\nu]\mathbf{m}_0$, because $Y_{00}(\theta, \phi) = 1/(4\pi)^{1/2}$. Since $B_{em}^{kn} = \delta_{ek}\delta_{mn}$, the final result is

$$\tilde{M}(\nu) = \frac{j_1(\nu)}{\nu} \mathbf{m}_0 + \sum_{\beta=1}^{2} g_{\beta} \sum_{k=0}^{\infty} \sum_{n=-k}^{k} (-i)^k c_{kn}^{\beta} A_{kk}^{\beta}(\nu) Y_{kn}(\theta_q, \phi_q). \qquad (107)$$

Using the properties of the coefficients $c_{kn}^\beta$ studied in Appendix A, we reformulate (107) in terms of the associated Legendre polynomials, the indices $k=2\nu$ and $n=2\mu$, and the coefficients $a_{\nu\mu}^\beta$:

$$\tilde{M}(\nu) = \frac{j_1(\nu)}{\nu} \mathbf{m}_0 + \sum_{\beta=1}^{2} g_{\beta} \sum_{\nu=1}^{\infty} \sum_{\mu=0}^{\nu} (-1)^{\nu} a_{\nu\mu}^{\beta} \rho_{\nu}^{\beta}(\nu) P_{2\nu}^{2\mu}(\cos \theta_{q}) \cos(2\mu\phi_{q}), \qquad (108)$$

where the radial function is given by

$$\rho_{\nu}^{\beta}(\nu) = - \frac{\nu j_{2\nu-1}(\nu) j_{2\nu}(i\kappa_{\beta}) - i\kappa_{\beta} j_{2\nu-1}(i\kappa_{\beta}) j_{2\nu}(\nu)}{\nu^2 + \kappa_{\beta}^2}. \qquad (109)$$

Using once again expression (70) and the well known series (71) for the Bessel functions of the first kind, we redefine the spherical Bessel functions of imaginary arguments [as they appear in (109)] as

$$\Upsilon_{\nu}(\tau) = j_{2\nu}(i\tau) = \frac{\pi^{1/2}}{2} \sum_{s=0}^{\infty} \frac{(-1)^{\nu}(\tau/2)^{2(s+\nu)}}{s!\Gamma(2\nu + s + 3/2)}, \qquad (110)$$

$$\aleph_{\nu}(\tau) = ij_{2\nu-1}(i\tau) = \frac{\pi^{1/2}}{2} \sum_{s=0}^{\infty} \frac{(-1)^{\nu}(\tau/2)^{2(s+\nu)-1}}{s!\Gamma(2\nu + s + 1/2)}, \qquad (111)$$

such that it becomes clear that $\tilde{M}(\nu)$ is a purely real-valued function. The radial function is then rewritten as

$$\rho_{\nu}^{\beta}(\nu) = - \frac{\nu j_{2\nu-1}(\nu) \Upsilon_{\nu}(\kappa_{\beta}) - \kappa_{\beta} \aleph_{\nu}(\kappa_{\beta}) j_{2\nu}(\nu)}{\nu^2 + \kappa_{\beta}^2}. \qquad (112)$$

By comparing this result with (102), we find the following pair of spherical Hankel transforms,

$$\rho_{\nu}^{\beta}(\nu) = \int_{0}^{1} \Upsilon_{\nu}(\kappa_{\beta}\xi) j_{\nu}(\nu\xi) \xi^2 d\xi. \qquad (113)$$

By comparing the result of (108) with (16) it becomes clear that the angular part is (due to the orthogonality relation of the spherical harmonics) shape invariant under Fourier transformation, while the spherical Hankel transform (113) of the radial function $\Upsilon_{\nu}(\kappa_{\beta}\xi)$ remains.

## Funding information

Funding for this research was provided by National Research Fund of Luxembourg (grant No. 15639149 to MPA and AM).

## References

Adams, M. P., Michels, A. & Kachkachi, H. (2022). *J. Appl. Cryst.* **55**, 1488–1499.

Battle, X., Moya, C., Escoda-Torroella, M., Iglesias, Ò., Fraile Rodríguez, A. & Labarta, A. (2022). *J. Magn. Magn. Mater.* **543**, 168594.

Bender, P., Honecker, D. & Fernández Barquín, L. (2019). *Appl. Phys. Lett.* **115**, 132406.

Berger, L., Labaye, Y., Tamine, M. & Coey, J. M. D. (2008). *Phys. Rev. B*, **77**, 104431.

Bersweiler, M., Bender, P., Vivas, L. G., Albino, M., Petrecca, M., Mühlbauer, S., Erokhin, S., Berkov, D., Sangregorio, C. & Michels, A. (2019). *Phys. Rev. B*, **100**, 144434.

Bertotti, G. (1998). *Hysteresis in Magnetism*. San Diego: Academic Press.

Brown, W. F. Jr (1963). *Micromagnetics*. New York: Interscience Publishers.

Dennis, C. L., Krycka, K. L., Borchers, J. A., Desautels, R. D., van Lierop, J., Huls, N. F., Jackson, A. J., Gruettner, C. & Ivkov, R. (2015). *Adv. Funct. Mater.* **25**, 4300–4311.

Disch, S., Wetterskog, E., Hermann, R. P., Wiedenmann, A., Vainio, U., Salazar-Alvarez, G., Bergström, L. & Brückel, T. (2012). *New J. Phys.* **14**, 013025.

Erokhin, S., Berkov, D. & Michels, A. (2015). *Phys. Rev. B*, **92**, 014427.

Garanin, D. A. & Kachkachi, H. (2003). *Phys. Rev. Lett.* **90**, 065504.

Garanin, D. A. & Kachkachi, H. (2009). *Phys. Rev. B*, **80**, 014420.

Gradmann, U. (1986). *J. Magn. Magn. Mater.* **54–57**, 733–736.

Gradshteyn, I. S. & Ryzhik, I. M. (2007). *Table of Integrals, Series, and Products*, edited by A. Jeffrey & D. Zwillinger. Burlington: Academic Press.

Grutter, A. J., Krycka, K. L., Tartakovskaya, E. V., Borchers, J. A., Reddy, K. S. M., Ortega, E., Ponce, A. & Stadler, B. J. H. (2017). *ACS Nano*, **11**, 8311–8319.
Günther, A., Bick, J.-P., Szary, P., Honecker, D., Dewhurst, C. D., Keiderling, U., Feoktystov, A. V., Tschöpe, A., Birringer, R. & Michels, A. (2014). *J. Appl. Cryst.* **47**, 992–998.

Hasz, K., Ijiri, Y., Krycka, K. L., Borchers, J. A., Booth, R. A., Oberdick, S. & Majetich, S. A. (2014). *Phys. Rev. B*, **90**, 180405.

Honecker, D., Bersweiler, M., Erokhin, S., Berkov, D., Chesnel, K., Venero, D. A., Qdemat, A., Disch, S., Jochum, J. K., Michels, A. & Bender, P. (2022). *Nanoscale Adv.* **4**, 1026–1059.

Honecker, D. & Michels, A. (2013). *Phys. Rev. B*, **87**, 224426.

Ijiri, Y., Krycka, K. L., Hunt-Isaak, I., Pan, H., Hsieh, J., Borchers, J. A., Rhyne, J. J., Oberdick, S. D., Abdelgawad, A. & Majetich, S. A. (2019). *Phys. Rev. B*, **99**, 094421.

Jackson, J. D. (1999). *Classical Electrodynamics*, 3rd ed. Hoboken: Wiley.

Kachkachi, H. (2007). *J. Magn. Magn. Mater.* **316**, 248–254.

Köhler, T., Feoktystov, A., Petracic, O., Kentzinger, E., Bhatnagar-Schöffmann, T., Feygenson, M., Nandakumaran, N., Landers, J., Wende, H., Cervellino, A., Rücker, U., Kovács, A., Dunin-Borkowski, R. E. & Brückel, T. (2021a). *Nanoscale*, **13**, 6965–6976.

Köhler, T., Feoktystov, A., Petracic, O., Nandakumaran, N., Cervellino, A. & Brückel, T. (2021b). *J. Appl. Cryst.* **54**, 1719–1729.

Krycka, K. L., Borchers, J. A., Booth, R. A., Ijiri, Y., Hasz, K., Rhyne, J. J. & Majetich, S. A. (2014). *Phys. Rev. Lett.* **113**, 147203.

Maurer, T., Gautrot, S., Ott, F., Chaboussant, G., Zighem, F., Cagnon, L. & Fruchart, O. (2014). *Phys. Rev. B*, **89**, 184423.

Metlov, K. L. & Michels, A. (2015). *Phys. Rev. B*, **91**, 054404.

Metlov, K. L. & Michels, A. (2016). *Sci. Rep.* **6**, 25055.

Mettus, D. & Michels, A. (2015). *J. Appl. Cryst.* **48**, 1437–1450.

Michels, A. (2021). *Magnetic Small-Angle Neutron Scattering: A Probe for Mesoscale Magnetism Analysis*. Oxford University Press.

Michels, A., Erokhin, S., Berkov, D. & Gorn, N. (2014). *J. Magn. Magn. Mater.* **350**, 55–68.

Michels, A., Mettus, D., Honecker, D. & Metlov, K. L. (2016). *Phys. Rev. B*, **94**, 054424.

Michels, A., Mettus, D., Titov, I., Malyeyev, A., Bersweiler, M., Bender, P., Peral, I., Birringer, R., Quan, Y., Hautle, P., Kohlbrecher, J., Honecker, D., Fernández, J. R., Barquín, L. F. & Metlov, K. L. (2019). *Phys. Rev. B*, **99**, 014416.

Mistonov, A. A., Dubitskiy, I. S., Shishkin, I. S., Grigoryeva, N. A., Heinemann, A., Sapoletova, N. A., Valkovskiy, G. A. & Grigoriev, S. V. (2019). *J. Magn. Magn. Mater.* **477**, 99–108.

Mühlbauer, S., Honecker, D., Périgo, A., Bergner, F., Disch, S., Heinemann, A., Erokhin, S., Berkov, D., Leighton, C., Eskildsen, M. R. & Michels, A. (2019). *Rev. Mod. Phys.* **91**, 015004.

Nedelkoski, Z., Kepaptsoglou, D., Lari, L., Wen, T., Booth, R. A., Oberdick, S. D., Galindo, P. L., Ramasse, Q. M., Evans, R. F. L., Majetich, S. & Lazarov, V. K. (2017). *Sci. Rep.* **7**, 45997.

Néel, L. (1954). *J. Phys. Radium*, **15**, 225–239.

Oberdick, S. D., Abdelgawad, A., Moya, C., Mesbahi-Vasey, S., Kepaptsoglou, D., Lazarov, V. K., Evans, R. F. L., Meilak, D., Skoropata, E., van Lierop, J., Hunt-Isaak, I., Pan, H., Ijiri, Y., Krycka, K. L., Borchers, J. A. & Majetich, S. A. (2018). *Sci. Rep.* **8**, 3425.

O’Handley, R. C. (2000). *Modern Magnetic Materials: Principles and Applications*. New York: Wiley.

Olver, F. W. J., Lozier, D. W., Boisvert, R. F. & Clark, C. W. (2010). *NIST Handbook of Mathematical Functions*. Cambridge University Press.

Pathak, S. A. & Hertel, R. (2021). *Phys. Rev. B*, **103**, 104414.

Riley, K. F., Hobson, M. P. & Bence, S. J. (2006). *Mathematical Methods for Physics and Engineering*. Cambridge University Press.

Vivas, L. G., Yanes, R., Berkov, D., Erokhin, S., Bersweiler, M., Honecker, D., Bender, P. & Michels, A. (2020). *Phys. Rev. Lett.* **125**, 117201.

Vivas, L. G., Yanes, R. & Michels, A. (2017). *Sci. Rep.* **7**, 13060.

Weber, H. J. & Arfken, G. B. (2003). *Essential Mathematical Methods for Physicists*. Amsterdam: Academic Press.

Wetterskog, E., Tai, C. W., Grins, J., Bergström, L. & Salazar-Alvarez, G. (2013). *ACS Nano*, **7**, 7132–7144.

Zákutná, D., Nižňanský, D., Barnsley, L. C., Babcock, E., Salhi, Z., Feoktystov, A., Honecker, D. & Disch, S. (2020). *Phys. Rev. X*, **10**, 031019.

Zaporozhets, V. D., Oba, Y., Michels, A. & Metlov, K. L. (2022). *J. Appl. Cryst.* **55**, 592–600.