# Fourier Transform of the Lippmann-Schwinger Equation:
Solving Vectorial Electromagnetic Scattering by
Arbitrary Shapes

Frédéric Gruy, Victor Rabiet, Mathias Perrin

► To cite this version:

Frédéric Gruy, Victor Rabiet, Mathias Perrin. Fourier Transform of the Lippmann-Schwinger Equation: Solving Vectorial Electromagnetic Scattering by Arbitrary Shapes. Mathematics, 2023, 11 (22), pp.4691. 10.3390/math11224691. emse-04355710

HAL Id: emse-04355710

https://hal-emse.ccsd.cnrs.fr/emse-04355710

Submitted on 20 Dec 2023

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
Article

# Fourier Transform of the Lippmann-Schwinger Equation: Solving Vectorial Electromagnetic Scattering by Arbitrary Shapes

Frederic Gruy<fnref n="1" />, Victor Rabiet <fnref n="1" />,<fnref n="2" /> and Mathias Perrin <fnref n="2" />,*

**Abstract:** In Electromagnetics, the field scattered by an ensemble of particles—of arbitrary size, shape, and material—can be obtained by solving the Lippmann-Schwinger equation. This singular vectorial integral equation is generally formulated in the direct space $\mathbb{R}^n$ (typically $n = 2$ or $n = 3$). In the article, we rigorously computed the Fourier transform of the vectorial Lippmann-Schwinger equation in the space of tempered distributions, $S'(\mathbb{R}^3)$, splitting it in a singular and a regular contribution. One eventually obtains a simple equation for the scattered field in the Fourier space. This permits to draw an explicit link between the shape of the scatterer and the field through the Fourier Transform of the body indicator function. We compare our results with accurate calculations based on the T-matrix method and find a good agreement.

**Keywords:** electromagnetic scattering; integral equation; singular integral; Fourier Transform

**MSC:** 46F12; 42A38; 31B20; 65R20

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

### check for updates

**Citation:** Gruy, F.; Rabiet, V.; Perrin, M. Fourier Transform of the Lippmann-Schwinger Equation: Solving Vectorial Electromagnetic Scattering by Arbitrary Shapes.

*Mathematics* **2023**, *11*, 4691. https://doi.org/10.3390/math11224691

Academic Editor: Nikolaos L. Tsitsas

Received: 28 September 2023

Revised: 13 November 2023

Accepted: 14 November 2023

Published: 18 November 2023

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

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

## 1. Introduction

The first use of Fourier transform (FT) to describe light scattering might be traced back to the early theories of Kirchhoff or Rayleigh and Sommerfeld [1]. The solution of Helmholtz equation that describes the field diffracted by an aperture in a screen was indeed related to the FT of its shape. Nowadays, solvers based on the Fourier expansion of electromagnetic fields are very efficient to model 1D or 2D periodic gratings [2,3], using the RCWA method [4,5], in close connection with experiments [6–8]. For 3D problems, a few results based on the FT of Maxwell equations are available. The calculation of the free space Green tensor [9–12], or a far-field description of the scattered field for a smooth particle [13], are presented. Using FT, a link between the far field scattered intensity and some geometrical properties of the scatterer (its Gaussian curvature) has been established, e.g., for the Porod’s law [14,15]. In other works, the polarization of electromagnetic waves in complex systems with scatterers was taken into account for an ensemble of point dipoles, e.g., solving the Bethe–Salpeter Equation in Fourier space [16]. However, to compute the vectorial scattered fields in arbitrary 3D geometries, the most popular methods [17], such as finite elements method (FEM) [18] or T-matrix computation [19], e.g., based on an expansion on vectorial spherical harmonics (VSH) [20–22], have been developed in real space.

These real space methods are, nevertheless, limited. For example, the memory requirement for FEM is huge when large spatial domains are to be modeled. The convergence of the T-matrix method, computed using the VSH ([17] see 3.1), is still debated [23] when close non-spherical scatterers—whose circumscribing spheres intersect—are studied. It is therefore still interesting to propose new analytical methods, which is what we shall do in this article.

<footnote n="1">¹ Ecole Nationale Supérieure des Mines de St. Etienne, Centre SPIN, F-42100 Saint-Etienne, France</footnote>

<footnote n="2">² Univ. Bordeaux, CNRS, LOMA, UMR 5798, F-33400 Talence, France</footnote>

<footnote>* Correspondence: mathias.perrin@u-bordeaux.fr</footnote>
Budko and co-authors [24,25] made the first move towards a Fourier Transform of the Lippman Schwinger Equation (LSE) for 3D geometries. This work was limited, however, to a scatterer described by Hölder continuous functions, and dealt with the singular part of the LSE only. It shall be interesting to go further and make the FT of the whole LSE. This would permit to (i) relate the field scattered by objects to the Fourier transform of the indicator function (The indicator function—whose value is 1 inside the scatterer(s), 0 outside—is known to be a convenient descriptor of the geometrical properties of the scatterer [14]); to factorize the equation in the Fourier space, instead of dealing with a convolution in real space, and (ii) to solve for the scattered field in fast and convenient way, using tridimensional Fourier transforms [26].

In Section 2, we conduct a rapid review of the regularization procedure that permits us to unambiguously define the Integral equation, see [27,28], and split the kernel in a regular and a singular contribution. Following that, we will detail the calculation of the FT of each part of the equation, see Sections 3 and 4, in the general case—as long as the scatterer is non-magnetic. For the singular part, see Section 3, we rely on [29] to explicitly express the FT. For the regular part, see Section 4, our approach is based on the work by Grafakos and Teschl [30] regarding the FT of radial (generalized) functions.

We shall gather the results in Section 5 to give the Fourier transform of electromagnetic LSE in $\mathbb{R}^3$. Eventually, we give a numerical example, solving the general equation in a few typical cases.

Note that the lengthy proofs and lemma have been gathered in [31].

## 2. Backgrounds and Purpose

*The LS Equation for Electromagnetic Scattering*

Stemming directly from the Maxwell equations, the LS equation that describes the field scattered by an arbitrary object in a uniform background is [24,28], in $\mathbb{R}^3$,

$$ \mathbf{E}(\mathbf{x}, \omega) = \mathbf{E}_{\text{inc}}(\mathbf{x}, \omega) + [k_b^2 + \text{grad div}] \int_{\mathbb{R}^3} G(\|\mathbf{x} - \mathbf{x}'\|) \chi(\mathbf{x}', \omega) \mathbf{E}(\mathbf{x}', \omega) \, d\mathbf{x}', \quad (1) $$

where

* $\mathbf{E}$ and $\mathbf{E}_{\text{inc}}$ are the total and incident electric field;

* $k_b$ is the wavenumber in the background medium, assumed to be real;

* $\chi$ is the relative difference between scatterer and background complex permittivity: $\chi(\mathbf{x}) := [\varepsilon(\mathbf{x}) - \varepsilon_b]/\varepsilon_b$. The scatterer medium is possibly dispersive—$\chi$ depends on $\omega$, dissipative—$\chi$ is complex—, and anisotropic—$\chi$ is a $3 \times 3$ tensor. In the case of an isotropic medium, $\chi$ boils down to a scalar. If $\varepsilon(\mathbf{x})$ is constant inside the scatterer, $\chi$ is proportional to the indicator function, and will be denoted like this in the following.

* $G$ is the Green function of Helmholtz equation in vacuum (with $l := \mathbf{x} - \mathbf{x}'$ and $l := \|\mathbf{1}\|$).

Note that if one chooses the $\exp(-i\omega t)$ convention for harmonic time evolution (as we do), the outgoing Sommerfeld condition [32] imposes:

$$ G(l) := \frac{e^{+ik_b l}}{4\pi l}, \quad (2) $$

so that far enough from the scatterer, the wave behaves as an outgoing spherical wave that carries energy outwards.
As shown in textbooks ([28], p. 36), ([33], chap. 15), the **grad** div operator in Equation (1) cannot be put straightforwardly under the integral sign. However, upon a careful procedure, Equation (1) can be recast as the singular integral equation

$$ \mathbf{E}_{\text{inc}}(\mathbf{x}, \omega) = \left[ \mathbb{I} + \frac{\chi}{3} \right] \mathbf{E}(\mathbf{x}, \omega) - \lim_{\epsilon \to 0} \int_{\mathbf{x}' \in \mathbb{R}^3 \setminus B(\mathbf{x}, \epsilon)} \mathbb{G}_0(\mathbf{x} - \mathbf{x}') \chi(\mathbf{x}', \omega) \mathbf{E}(\mathbf{x}', \omega) \, d\mathbf{x}' \\ - \int_{\mathbf{x}' \in \mathbb{R}^3} \mathbb{G}_1(\mathbf{x} - \mathbf{x}') \chi(\mathbf{x}', \omega) \mathbf{E}(\mathbf{x}', \omega) \, d\mathbf{x}', \qquad (3) $$

where $\mathbb{I}$ is the unit tensor, $\mathbb{G}_0$ and $\mathbb{G}_1$ are, respectively, the “singular” part and the “regular” part of the Green tensors (with $\mathbb{I} := \mathbf{x} - \mathbf{x}'$ and $l := \|\mathbb{I}\|$):

$$ \mathbb{G}_0(\mathbb{I}) := \frac{1}{4\pi l^3} (3\mathbb{Q} - \mathbb{I}), \qquad (4) \\ \mathbb{G}_1(\mathbb{I}) := G_{1I}(l)\mathbb{I} + G_{1Q}(l)\mathbb{Q}, $$

where

$$ \mathbb{Q} := \frac{\mathbb{I}(\mathbb{I})}{l^2} = \frac{1}{l^2} \begin{pmatrix} l_x l_x & l_x l_y & l_x l_z \\ l_y l_x & l_y l_y & l_y l_z \\ l_z l_x & l_z l_y & l_z l_z \end{pmatrix}, \qquad (5) $$

and

$$ G_{1I}(l) := \frac{e^{ik_b l}}{4\pi l^3} (-1 + ik_b l - (ik_b l)^2) + \frac{1}{4\pi l^3}, \qquad (6) $$

$$ G_{1Q}(l) := \frac{e^{ik_b l}}{4\pi l^3} (3 - 3ik_b l + (ik_b l)^2) - \frac{3}{4\pi l^3}. \qquad (7) $$

Note that in the quasistatic limit, i.e., when $1/k_b$ is much larger than the size of the scatterer, the term with $\mathbb{G}_1$ cancels this out and only the principal value that contains $\mathbb{G}_0$ contributes.

Equation (3) defines a linear operator whose properties have been studied both in 2D [25] and 3D [24].

Our purpose is to give an explicit formula for the Fourier transform of Equation (3), with the key point being the computation of the Fourier transform of Green tensors $\mathbb{G}_0$ and $\mathbb{G}_1$. We will proceed in two steps. $\mathbb{G}_0$ will be Fourier transformed in a more general way than what has previously been carried out [24], and $\mathbb{G}_1$ will be handled using generalized functions theory.

### 3. Fourier Transform of the Singular Part

The definitions and properties of the Fourier Transform are given in Appendix A. There, we introduce a general notation for the FT, that depends on two parameters $a$ and $b$, see Equation (A1). In the article, all the calculations are carried using $(a, b) = (1, -1)$. The results in the general case are then given without demonstration so that they can be easily adapted to any definition of the FT.

#### 3.1. General Expression as a Convolution

Let $f$ be an integrable function on the sphere $S^{n-1}$ with mean 0. Following that, the FT can be defined, possibly as a generalized function.

We define accordingly

$$ \langle W_f, \varphi \rangle := \lim_{\epsilon \to 0} \int_{\mathbb{R}^n \setminus B(0, \epsilon)} \frac{f(x/|x|)}{|x|^n} \varphi(x) \, dx, \qquad (8) $$

for $\varphi \in S(\mathbb{R}^d)$. It has been shown that $W_f$ is a tempered distribution (see e.g., [34], p. 334), provided that the hypothesis “$f$ has a zero mean over the sphere $S$” is true (see
also [24]). Eventually, the Fourier transform can be found (cf. [34], p. 336). It is the function (i.e., finite) given by:

$$ \widehat{W}_f(\xi) = \int_{\mathbb{S}^{n-1}} f(\theta) \left( \log \frac{1}{|\xi \cdot \theta|} - \frac{i\pi}{2} \operatorname{sgn}(\xi \cdot \theta) \right) d\theta, \quad (9) $$

where $\mathbb{S}^{n-1}$ is the unit sphere on $\mathbb{R}^n$, and $\theta \in \mathbb{S}^{n-1}$ is the direction of the vector $x/|x|$ for $x \in \mathbb{R}^n$.

To carry on with the Fourier transform of the singular part, we defined the continuous operator $H_f$ by

$$ H_f(u) := W_{f,1} \star u + W_{f,2} \star u, \quad (10) $$

where $\star$ represents the convolution, with $u \in L^2$, and

$$ \langle W_{f,1}, u \rangle := \lim_{\epsilon \to 0} \int_{B(0,1) \setminus B(0,\epsilon)} \frac{f(x/|x|)}{|x|^n} u(x) dx \quad \text{and} \quad W_{f,2} := 1_{B(0,1)^c} \frac{f(x/|x|)}{|x|^n}, $$

where $B(0,1)^c$ is the complementary of the ball of radius 1 in $\mathbb{R}^n$.

Besides, we will assume the following hypothesis:

**Hypothesis 1.** $\widehat{W}_f$ is bounded.

**Remark 1.** As we will see, that will be always the case in our framework.

**Theorem 1.** $H_f$ is a continuous operator over the Hilbert space $L^2(\mathbb{R}^n)$ with a Fourier transform verifying:

$$ \widehat{H}_f(u) = \widehat{W}_f \hat{u} \in L^2. \quad (11) $$

**Proof.** The operator $H_f$ is well-defined since, on one hand, $W_{f,1}$ is a distribution with compact support, $u$ can be assimilated to a tempered distribution ($L^2 \subset S'$) and, on the other hand $W_{f,2}$ and $u$ are in $L^2$.

Additionally, the Fourier transform of $H_f$ is well-defined and

$$ \begin{aligned} \widehat{H}_f(u) &= \widehat{W}_{f,1} \hat{u} + \widehat{W}_{f,2} \hat{u} \\ &= \widehat{W}_f \hat{u} \in L^2, \end{aligned} $$

since $\widehat{W}_f$ is supposed to be bounded. Then $\bar{\mathcal{F}}(\widehat{H}_f) \in L^2$, where $\bar{\mathcal{F}}$ denotes the inverse Fourier Transform.

The linear behaviour of $H_f$ is clear and, for the continuity, we have

$$ \|H_f(u)\|_2 = \|\widehat{H}_f(u)\|_2 = \|\widehat{W}_f \hat{u}\|_2 \le C \|\hat{u}\|_2 = C \|u\|_2. $$

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

This theorem shows that one can factorize the singular part of the LS equation Equation (3), in the Fourier Domain. Such a factorization was obtained by Budko [24] in the case of a Hölder continuous indicator function. Theorem 1 permits to generalize it to step index profiles (with, e.g., an indicator that is a Heaviside function), for which the function $x' \mapsto \chi(x', \omega) \cdot E(x', \omega)$ is $L^2$ but not Hölder.

We now have to compute $\widehat{W}_f$.

## 3.2. Explicit Calculation

In agreement with the definition of $G_0$, let us write each component of $f$ as an homogeneous polynomial $P(\theta) = c \prod_{k=1}^n \theta_k^{\alpha_k}$, where $\alpha_k$ are, respectively, even or odd numbers.
First, we can simply say that the part with the log function, in Equation (9), is null if $f$ is an odd function, and the part with the sgn function is null if $f$ is an even function.

$P(\theta)$ can be seen as a component of the polyadic tensor of order $s = \sum_{k=1}^{n} \alpha_k$, whose FT can be deduced from [29].

The obtained results for low order tensors (up to 2) are gathered in the following tables, where $\Psi$ denotes the digamma function. Detail on the calculation of the integrals is provided in Appendix B.2, whereas the general demonstration can be found in [29].

In passing, we note that for polyadic tensors $W_f$, whose components are given in Table 1, the Fourier Transformed tensor, $\hat{W}_f$ are bounded operators, as seen from their components (see Tables 1 and 2). This corresponds to our working hypothesis.

<table><caption>**Table 1.** Singular integrals Equation (9) for even and odd tensors.</caption><thead><tr><th>Tensor Components</th><th>Components of W<sub>f</sub>(ξ)</th></tr></thead><tbody><tr><td>1</td><td>E'<sub>0</sub></td></tr><tr><td>θ<sub>l</sub></td><td>O<sub>l</sub> = O'<sub>1</sub>ξ<sub>l</sub></td></tr><tr><td>θ<sub>l</sub>θ<sub>m</sub></td><td>E<sub>lm</sub> = E'<sub>11</sub>ξ<sub>l</sub>ξ<sub>m</sub> + E'<sub>22</sub>(-ξ<sub>l</sub>ξ<sub>m</sub> + δ<sub>l,m</sub>)</td></tr></tbody></table>

<table><caption>**Table 2.** Coefficients for the singular integrals of the even and odd tensors.</caption><thead><tr><th>Tensor Order</th><th>Definition of Coefficients</th><th>Value</th></tr></thead><tbody><tr><td>0</td><td>E'<sub>0</sub> = -∫<sub>∂B(0,1)</sub> log(|θ'<sub>1</sub>|) dθ'</td><td>π<sup>n/2</sup> <sup>Ψ(n/2)-Ψ(1/2)</sup>&frasl;<sub>Γ(n/2)</sub></td></tr><tr><td>1</td><td>O'<sub>1</sub> = -iπ/∫<sub>θ'<sub>1</sub>>0</sub> θ'<sub>1</sub> dθ'</td><td>-(iπ/2)π<sup>n-1/2</sup> <sup>2</sup>&frasl;<sub>Γ((n+1)/2)</sub></td></tr><tr><td rowspan="2">2</td><td>E'<sub>11</sub> = -∫<sub>∂B(0,1)</sub> log(|θ'<sub>1</sub>|)θ'<sub>1</sub>θ'<sub>1</sub> dθ'</td><td>π<sup>n/2</sup> <sup>1</sup>&frasl;<sub>2</sub> <sup>Ψ((n+2)/2)-Ψ(3/2)</sup>&frasl;<sub>Γ((n+2)/2)</sub></td></tr><tr><td>E'<sub>22</sub> = -∫<sub>∂B(0,1)</sub> log(|θ'<sub>1</sub>|)θ'<sub>2</sub>θ'<sub>2</sub> dθ'</td><td>π<sup>n/2</sup> <sup>1</sup>&frasl;<sub>2</sub> <sup>Ψ((n+2)/2)-Ψ(1/2)</sup>&frasl;<sub>Γ((n+2)/2)</sub></td></tr></tbody></table>

### 3.2.1. The Tridimensional Case, $n = 3$

The case of $\mathbb{R}^3$ is of particular interest. Let us give some details, using the results:

$$ \int_{\mathbb{S}^2} \log |\theta_1| d\theta = -4\pi, \quad \int_{\mathbb{S}^2} \log |\theta_1|\theta_1^2 d\theta = -\frac{4}{9}\pi, \quad \int_{\mathbb{S}^2} \log |\theta_1|\theta_2^2 d\theta = -\frac{16}{9}\pi. $$

From Equation (9), one obtains that, if $f = \eta$ is a polynomial of degree 0,

$$ \hat{W}_{\eta}(\xi) = 4\pi \eta. \qquad (12) $$

Now, if $f((x_1, x_2, x_3)) = \eta x_i x_j$,

$$ \hat{W}_f(\xi) = -\frac{4\pi}{3}f(\xi) + \frac{16\pi}{9}\eta \delta_{i,j}. \qquad (13) $$

So, if $f(l) = \frac{1}{4\pi}(3\mathbb{Q} - 1)$,

$$ \hat{W}_f(\xi) = F(\hat{\Theta}) = \frac{1}{3}\mathbb{I} - \hat{\mathbb{Q}}, \qquad (14) $$

where $\hat{\Theta} := \frac{\xi}{\|\xi\|}$, and $\hat{\mathbb{Q}} := \frac{\xi^{(l)\xi}}{\|\xi\|^2}$. The notation $F(\hat{\Theta})$ is introduced in agreement with [24], to emphasize that, as $f$ is of zero mean on $\mathbb{S}^{n-1}$, the Fourier transform of the singular part only depends on $\frac{\xi}{\|\xi\|}$.
### 3.2.2. General Notation

Note that we can easily express the result in the general notation, as

$$ \hat{W}_f(\xi) = F(\tilde{\Theta}) = c_{a,b}^3 \left( \frac{1}{3} \mathbb{I} - \hat{\mathbb{Q}} \right), \qquad (15) $$

where $c_{a,b}^3$ is defined as

$$ c_{a,b}^3 = \left( \frac{|b|}{(2\pi)^{1-a}} \right)^{\frac{3}{2}}. \qquad (16) $$

## 4. Fourier Transform of the Regular Part

### 4.1. Radial Distributions and Grafakos Theorem

As we will see below, in order to compute the FT of the regular part of Equation (3), we shall compute FT of radial distributions, and use specific results on this topic, demonstrated in [30]. Let us first recall some general definitions and properties. Following that, we will limit ourselves to a space dimension $n = 3$.

#### 4.1.1. Radial Distributions

While $S(\mathbb{R}^n)$ stands for the space of Schwartz functions on $\mathbb{R}^n$, we set:

$$ S_{\text{rad}}(\mathbb{R}^n) = \{ \varphi \in S(\mathbb{R}^n) : \varphi = \varphi \circ A, \forall A \in \mathcal{O}(n) \} \qquad (17) $$

$$ S_{\text{rad}}(\mathbb{R}) = S_{\text{even}}(\mathbb{R}) = \{ \varphi \in S(\mathbb{R}) : \varphi(x) = \varphi(-x) \} \qquad (18) $$

where $\mathcal{O}(n)$ is the set of the orthogonal transformations of $\mathbb{R}^n$.

We define then the following functions:

$$ \left\{ \begin{array}{l} S(\mathbb{R}^n) \rightarrow S_{\text{rad}}(\mathbb{R}) \\ \varphi \mapsto \left( r \mapsto \varphi^o(r) := \frac{1}{\omega_{n-1}} \int_{S^{n-1}} \varphi(r\theta) \, d\theta \right) \end{array} \right. \qquad (19) $$

$$ \left\{ \begin{array}{l} S_{\text{rad}}(\mathbb{R}) \rightarrow S_{\text{rad}}(\mathbb{R}^n) \\ \varphi \mapsto (x \mapsto \varphi^O(x) := \varphi(|x|)) \end{array} \right. \qquad (20) $$

(where $S^{n-1}$ is the unit sphere on $\mathbb{R}^n$ and $\omega_{n-1}$ its surface area; with the convention $\omega_0 = 2$ and $\varphi^o(x) = \frac{1}{2}(\varphi(x) + \varphi(-x))$, for $\varphi \in S(\mathbb{R})$).

**Definition 1.** A distribution $u \in S'(\mathbb{R}^n)$ (with $S'(\mathbb{R}^n)$ the space of tempered distributions on $\mathbb{R}^n$) is called radial if for all $A \in \mathcal{O}(n)$,

$$ u = u \circ A, $$

that is,

$$ \langle u, \varphi \rangle = \langle u, \varphi \circ A \rangle, $$

for all $\varphi \in S(\mathbb{R}^n)$. The set of all radial tempered distributions is denoted by $S'_{\text{rad}}(\mathbb{R}^n)$.

**Proposition 1.** For $u \in S'_{\text{rad}}(\mathbb{R}^n)$ and $\varphi \in S(\mathbb{R}^n)$,

$$ \langle u, \varphi \rangle = \langle u, \varphi^{\text{rad}} \rangle, \qquad (21) $$

where $\varphi^{\text{rad}} := (\varphi^o)^\circ$ (i.e., $\varphi^{\text{rad}}(x) = \varphi^o(|x|)$).

Given a function $\varphi : \mathbb{R}^n \rightarrow \mathbb{R}$ and a fixed $x \in \mathbb{R}^n$, $\varphi^o(|x|)$ is then the mean over the sphere of radius $|x|$ of the function $\varphi$ (so, if $\varphi$ is already radial, $\varphi$ will be of the form $f(|x|)$ with $f : \mathbb{R} \rightarrow \mathbb{R}$, and we will have directly $\varphi^o = f$).
Let us define the space

$$ \mathcal{R}_n := r^{n-1}\mathcal{S}_{\text{rad}}(\mathbb{R}) = \{(r \mapsto \psi(r)r^{n-1}), \psi \in \mathcal{S}_{\text{rad}}(\mathbb{R})\}. \quad (22) $$

**Remark 2.** $\mathcal{R}_n$ is a subspace $\mathcal{S}(\mathbb{R})$ on which we can use the same topology; we denote its dual (set of the linear continuous functions defined over $\mathcal{R}_n$) by $\mathcal{R}_n'$.

We switched from $\mathcal{R}_n'$ to radials distributions of $\mathcal{S}'(\mathbb{R}^n)$ as follows:

* if $u$ is radial distribution, we define $u_{\diamond} \in \mathcal{R}_n'$ by

$$ \langle u_{\diamond}, \psi(r)r^{n-1} \rangle := \frac{2}{\omega_{n-1}} \langle u, \psi^{\circ} \rangle, \quad \psi \in \mathcal{S}_{\text{rad}}(\mathbb{R}) \quad (23) $$

* if $u_{\diamond} \in \mathcal{R}_n'$, we define a radial distribution $u$ by

$$ \langle u, \varphi \rangle := \frac{\omega_{n-1}}{2} \langle u_{\diamond}, \varphi^{\circ}(r) r^{n-1} \rangle, \quad \varphi \in \mathcal{S}(\mathbb{R}^n) \quad (24) $$

### 4.1.2. Grafakos–Teschl Theorem

**Theorem 2 ([30]).** Given $v_1$ in $\mathcal{S}'(\mathbb{R})$, we define a radial distribution $v_k$ on $\mathbb{R}^k$ ($k \in \mathbb{N}^*$) by

$$ \langle v_k, \varphi \rangle := \frac{\omega_{k-1}}{2} \langle v_1, \varphi^0(r) r^{k-1} \rangle, \quad \varphi \in \mathcal{S}_{\text{rad}}(\mathbb{R}^k) \quad (25) $$

(if $\varphi \in \mathcal{S}_{\text{rad}}(\mathbb{R}^n)$, then $\varphi(x) = \varphi^0(|x|)$).

Let $u^k = \mathcal{F}_k^{a,b}(v_k)$. We have then

$$ -\frac{(2\pi)^a}{|b|r} \frac{d}{dr} u_{\diamond}^n = u_{\diamond}^{n+2}. \quad (26) $$

In the following, we will limit ourselves to $n = 3$.

Besides, instead of using directly the Grafakos–Teschl theorem, we will use the following corollary:

**Corollary 1.** Given $v_1$ in $\mathcal{S}'(\mathbb{R})$, we define the radial distribution $v_3$ on $\mathbb{R}^3$ by

$$ \langle v_3, \varphi \rangle := \frac{\omega_2}{2} \langle v_1, \varphi^0(r) r^2 \rangle, \quad \varphi \in \mathcal{S}_{\text{rad}}(\mathbb{R}^3) \quad (27) $$

(if $\varphi \in \mathcal{S}_{\text{rad}}(\mathbb{R}^n)$, then $\varphi(x) = \varphi^0(|x|)$). We have then, for all $\varphi \in \mathcal{S}_{\text{rad}}(\mathbb{R}^3)$

$$ \langle \mathcal{F}_3^{a,b}(v_3), \varphi \rangle = -\frac{(2\pi)^{a+1}}{|b|} \langle r \frac{d}{dr}(\mathcal{F}_1^{a,b}(v_1)), \varphi^0(r) \rangle. \quad (28) $$

**Proof.** With $u^3 := \mathcal{F}_3^{a,b}(v_3)$ and $u^1 := \mathcal{F}_1^{a,b}(v_3)$ we have

$$ \langle u^3, \varphi \rangle = \frac{\omega_2}{2} \langle u_{\diamond}^3, \varphi^0(r) r^2 \rangle \quad \text{cf. (24)} $$

$$ = \frac{\omega_2}{2} \left\langle -\frac{(2\pi)^a}{|b|r} \frac{d}{dr} u^1, \varphi^0(r) r^2 \right\rangle \quad \text{cf. (26)} $$

$$ = -\frac{(2\pi)^{a+1}}{|b|} \langle r \frac{d}{dr}(u^1), \varphi^0(r) \rangle, $$

since $\omega_2 = 4\pi$. $\square$
## 4.2. Fourier Transform of the Regular Part of $G_{11}(l)\mathbb{I}$

### Procedure for Computing the FT in $\mathbb{R}^3$

We aim at computing the FT of a radial distribution, say $v_3$. We will proceed using Corollary (1), in the following way:

1. We need to find a one-dimensional distribution, $v_1$, which verifies

$$ \langle v_3, \varphi \rangle = \frac{\omega_2}{2} \langle v_1, \varphi^0(r) r^2 \rangle. \qquad (29) $$

($\omega_2 = 4\pi$ being the 1-sphere surface). Note that finding $v_1$ may be, generally speaking, a difficult task.

2. One computes the one-dimensional Fourier Transform $u_1 := \hat{v}_1$.

3. The FT we seek, $u_3 := \hat{v}_3$, is then defined by

$$ \langle u_3, \varphi \rangle = - \left( \frac{\omega_2}{2} \right)^2 \langle r \frac{du_1}{dr}, \varphi^0(r) \rangle. \qquad (30) $$

4. Note that if $\frac{du_1}{dr} = f(r)$ is not a distribution, but a regular function, one directly writes

$$ u_3 = -\pi \frac{f(\|l\|)}{\|l\|}. \qquad (31) $$

In practice, we want to compute the 3D Fourier transform of

$$ l \mapsto a_0 \frac{(e^{ik_b l} - 1)}{4\pi l^3} + \frac{e^{ik_b l}}{4\pi l^3} (a_1 l + a_2 l^2). \qquad (32) $$

And we will later impose, to correspond to the Maxwell equation, that

$$ a_0 = -1 $$

$$ a_1 = ik_b $$

$$ a_2 = k_b^2 $$

The idea is thus to reduce the problem to a one dimensional Fourier transform computation. To proceed, we will use Theorem 2. However, this is not sufficient. Indeed, in this theorem, a one-dimensional distribution (denoted by $v_1$) is given and is used to compute the associated (Equation (27)) multidimensional distribution $v_3$. In our work, on the contrary, we start from a multidimensional distribution $v_3$, and $v_1$ is not known a priori. To overcome this problem, we will first define (in a heuristic way) a one-dimensional distribution and use Theorem 2 to obtain an associated multidimensional distribution $T$ (using Equation (27)) whose Fourier Transform is established by this very theorem. Finally, to conclude, it remains necessary to prove that the Fourier transform of $T$ is the same as the Fourier transform of $v_3$. The proof relies on a recent result [35]), and we will not elaborate here.

Rewriting Equation (32), in the sense of generalized functions of $S'(\mathbb{R})$, we define $v_1$ as:

$$ v_1(r) = a_0 \left[ (\cos(k_b r) - 1) T_{\frac{1}{|r|^3}} + i \sin(k_b r) T_{\frac{1}{r^3}} \right] + a_1 \left[ \cos(k_b r) T_{\frac{1}{r^2}} + i \sin(k_b r) T_{\frac{\operatorname{sgn}(r)}{r^2}} \right] \\ + a_2 \left[ \cos(k_b r) T_{\frac{1}{|r|}} + i \sin(k_b r) T_{\frac{1}{r}} \right]. \qquad (33) $$

where $T_{\frac{1}{r^2}}$, $T_{\frac{\operatorname{sgn}(r)}{r^2}}$, and $T_{\frac{1}{|r|^3}}$ are the distributions associated to $r \mapsto \frac{1}{r^2}$, $r \mapsto \frac{\operatorname{sgn}(r)}{r^2}$ and $r \mapsto \frac{1}{|r|^3}$, respectively.

We start by computing the one-dimensional Fourier transform of $v_1$.
Using the previous results on the Fourier Transform of distributions ([31], see Sec. I of Supplementary Information), ref. [36] we found, without difficulty:

$$ \hat{v}_1 = \sum_{s \in \{-1,0,1\}} \left[ \gamma + \ln|k - sk_b| + i \frac{s\pi}{2} \operatorname{sgn}(k - sk_b) \right] \times \left[ A_2(k - sk_b)^2 + A_1(k - sk_b) + A_0 \right], \quad (34) $$

with

$$ \begin{aligned} A_0 &= -a_2 s^2, \\ A_1 &= ia_1 s, \\ A_2 &= a_0(3s^2 - 2)/2. \end{aligned} \quad (35) $$

Note that the detailed and rigorous computation is given in an external document for the specific and most important case of $n = 3$, and Maxwell equations. See, in particular, Ref. [31] (Sec. (IV. B) of the main text and Sec. (IV A.) of the Supplementary Information). To compute the FT of the associated 3D distribution $v_3$, according to Corollary 1 (with $a = 1, b = -1$), we will have to compute the derivative of $\hat{v}_1$, in the sense of generalized functions.

Using that

$$ \frac{d}{dk} \operatorname{sgn}(k - sk_b) = 2\delta(k - sk_b), \quad (36) $$

and the equalities

$$ (k - sk_b)\delta(k - sk_b) = 0, \quad (k - sk_b)^2\delta(k - sk_b) = 0, $$

we find that

$$ \begin{aligned} \frac{d\hat{v}_1}{dk} = & \sum_{s \in \{-1,0,1\}} A_2 (1 + 2\gamma)(k - sk_b) + A_1(1 + \gamma) + \frac{A_0}{k - sk_b} + \\ & A_1 \ln|k - sk_b| + 2A_2(k - sk_b) \ln|k - sk_b| + \\ & i \frac{\pi}{2} s [2A_0\delta(k - sk_b) + (A_1 + 2A_2(k - sk_b))\operatorname{sgn}(k - sk_b)]. \end{aligned} \quad (37) $$

Computing the sum explicitly, one obtains:

$$ \begin{aligned} \frac{d\hat{v}_1}{dk} = & -a_2 \left( T_{\frac{1}{k-k_b}} + T_{\frac{1}{k+k_b}} \right) + a_0 k \ln \frac{|k^2 - k_b^2|}{k^2} + (ia_1 - a_0 k_b) \ln \frac{|k - k_b|}{|k + k_b|} \\ & +i \frac{\pi}{2} \left[ a_0[k \operatorname{sgn}(k - k_b) - k \operatorname{sgn}(k + k_b)] + (ia_1 - a_0 k_b)[\operatorname{sgn}(k - k_b) + \operatorname{sgn}(k + k_b)] \right. \\ & \left. \qquad -2a_2(\delta(k - k_b) - \delta(k + k_b)) \right]. \end{aligned} \quad (38) $$

Then, using Corollary 1, with $a = 1, b = -1, u_3 = \hat{v}_3$ and $u_1 = \hat{v}_1$, we seek $\hat{v}_3$ such as

$$ \langle \hat{v}_3, \varphi \rangle = -\frac{\omega_2^2}{4} \left\langle k \frac{d}{dk} \hat{v}_1, \varphi^o(k) \right\rangle. \quad (39) $$

By linearity, $\hat{v}_3$ is the sum of the contribution of each term of $k \frac{d}{dk} \hat{v}_1$, calculated using Equation (39):
* The contribution of $-a_2k(T_{\frac{1}{k-k_b}} + T_{\frac{1}{k+k_b}})$ to $\hat{v}_3$ has been calculated, (see [31], Sec. IV. A. 3. of the Supplementary Information), and is

$$ \frac{\omega_2}{2} a_2 k^{-1} \left[ v.p. \left( \frac{1}{k-k_b} \right) - \frac{1}{k+k_b} \right], $$

* Rewriting $a_0k^2 \ln \frac{|k^2-k_b^2|}{k^2} = a_0k^2(\ln|k-k_b| + \ln|k+k_b| - 2\ln|k|)$, and using ([31], 1. and 2. page 14 of the main text), the contribution of this term to $\hat{v}_3$ is

$$ -\frac{\omega_2}{2} a_0 \ln \frac{|k^2-k_b^2|}{k^2}, $$

* Similarly, the term $(ia_1 - a_0k_b)k \ln \frac{|k-k_b|}{|k+k_b|}$ gives a contribution

$$ -\frac{\omega_2}{2} (ia_1 - a_0k_b)G(k) \text{ to } \hat{v}_3, $$

* the term $i\frac{\pi}{2}a_0k[k \operatorname{sgn}(k-k_b) - k \operatorname{sgn}(k+k_b)]$ gives a contribution, see ([31], 1. page 15 of the main text):

$$ i\pi a_0 \frac{\omega_2}{2} \mathbb{1}_{B(0,k_b)}(k), $$

* a similar calculation proves that the term $i\frac{\pi}{2}(ia_1 - a_0k_b)k[\operatorname{sgn}(k-k_b) + \operatorname{sgn}(k+k_b)]$ gives a contribution

$$ -i\pi \frac{\omega_2}{2} (ia_1 - a_0k_b)F(k), $$

* the term $-i\pi a_2k(\delta(k-k_b) - \delta(k+k_b))$ gives a contribution, (see [31], 2. page 15 of the main text):

$$ i\pi \frac{a_2}{k_b} \frac{\omega_2}{2} \delta(|k| - k_b), $$

where

$$ F(k) = H(|k| - k_b) / |k|, \text{ and } G(k) = k^{-1} \ln \left| \frac{k - k_b}{k + k_b} \right|, $$

which are both even functions of $k$.

Using the fact that

$$ \hat{v}_3 = \omega_2 \hat{G}_{1I}, $$

we obtain:

**Proposition 2.** The Fourier Transform of the regular part:

$$ l \mapsto \frac{e^{ik_b l}}{4\pi l^3} (a_0 + a_1 l + a_2 l^2) - \frac{a_0}{4\pi l^3} $$

is (with $k = \sqrt{k_x^2 + k_y^2 + k_z^2}$)

$$ k \mapsto a_0 P(k) + \frac{a_2}{k_b^2} A(k) - \frac{ia_1 - a_0 k_b}{2} \frac{M(k)}{k}, \qquad (40) $$
with the definitions:

$$A(k) = \frac{k_b^2}{k^2 - k_b^2} + i \frac{\pi}{2} k_b \delta(|k| - k_b), \qquad (41)$$

$$P(k) = -\frac{1}{2} \ln \left| \frac{k^2 - k_b^2}{k^2} \right| + i \frac{\pi}{2} \mathbb{1}_{B(0, k_b)}, \qquad (42)$$

$$M(k) = k(G(k) + i\pi F(k)) = \ln\left|\frac{k - k_b}{k + k_b}\right| + i\pi H(|k| - k_b). \qquad (43)$$

### 4.3. Fourier Transform of the Regular Part of $G_{1Q}(l)Q$

We now want to Fourier Transform the tensorial valued function

$$l \mapsto \left[ \frac{e^{ik_b l}}{4\pi l^3} (c_0 + c_1 l + c_2 l^2) - \frac{c_0}{4\pi l^3} \right] Q. \qquad (44)$$

And we will later impose

$$\begin{aligned} c_0 &= 3 \\ c_1 &= -3ik_b \\ c_2 &= -k_b^2 \end{aligned}$$

Let us define $l \mapsto G_{1Q}(l)$ with

$$G_{1Q}(l) = \left[ \frac{e^{ik_b l}}{4\pi l^3} (c_0 + c_1 l + c_2 l^2) - \frac{c_0}{4\pi l^3} \right]. \qquad (45)$$

Contrary to the previous case, the components of $G_{1Q}(l)Q$ do not possess radial symmetry, since they have the form $\frac{l_i l_j}{l^2} G_{1Q}(l) \ (\in L_{\text{loc}}^1(\mathbb{R}^3))$: we cannot directly use Corollary 1. Therefore, in order to proceed, we will first build a radial distribution $T$ such that

$$l_i l_j T = \frac{l_i l_j}{l^2} G_{1Q}(l).$$

We will then have (denoting by $D$ the differentiation in the sense of the generalized functions)

$$\mathcal{F}\left(\frac{l_i l_j}{l^2} G_{1Q}(l)\right) = \mathcal{F}(l_i l_j T) = -D_{k_i} D_{k_j} \mathcal{F}(T). \qquad (46)$$

Heuristically (since the following $\frac{1}{l^2} G_{1Q}(l)$ function is not definite at 0), we have to compute the 3D Fourier transform of $l \mapsto v_1(l)$ defined with,

$$v_1(l) = \left[ \frac{e^{ik_b l} - 1}{l^5} c_0 + \frac{e^{ik_b l}}{l^5} (c_1 l + c_2 l^2) \right],$$

and to use a unidimensional Fourier Transform computation (Corollary 1).

#### 4.3.1. Fourier Transform of $v_1$

Expressing $v_1$ as a generalized function, we obtain:

$$v_1(r) = c_0 \left[ \left( \cos(k_b r) - 1 \right) T_{\frac{1}{|r|^5}} + i \sin(k_b r) T_{\frac{1}{r^4}} \right] + c_1 \left[ \cos(k_b r) T_{\frac{1}{r^4}} + i \sin(k_b r) T_{\frac{\operatorname{sgn}(r)}{r^4}} \right] \\ + c_2 \left[ \cos(k_b r) T_{\frac{1}{|r|^3}} + i \sin(k_b r) T_{\frac{1}{r^3}} \right]. \qquad (47)$$
As in the previous case, we compute the one-dimensional Fourier Transform of $v_1$ in the sense of generalized functions, and obtain:

$$ \hat{v}_1 = \frac{-1}{2} \sum_{s \in \{-1,0,1\}} \left[ \gamma + \ln|k - sk_b| + i \frac{s\pi}{2} \operatorname{sgn}(k - sk_b) \right] \times \left[ \frac{C_2}{6}(k - sk_b)^4 + \frac{C_1}{3}(k - sk_b)^3 + C_0(k - sk_b)^2 \right], \quad (48) $$

where the $C_{0,1,2}$ are given by

$$ \begin{aligned} C_0 &= -c_2 s^2, \\ C_1 &= ic_1 s, \\ C_2 &= c_0 (3s^2 - 2)/2 \end{aligned} \quad (49) $$

Note that the detailed and rigorous computation is given in an external document for the specific case of $n = 3$, and Maxwell equations. See, in particular, [31], Sec. (IV. C) of main text and Sec. (IV B.) of Supplementary Information.

Following that, we need to compute the derivative of these expressions (in the sense of the generalized functions) using the Lemma (II.12), see ([31], p. 13 of Supplementary Information).

We eventually obtain:

$$ \frac{d\hat{v}_1}{dk} = -\frac{1}{2} \sum_{s \in \{-1,0,1\}} C_2(\frac{2}{3}\gamma + \frac{1}{6})(k - sk_b)^3 + C_1(\gamma + \frac{1}{3})(k - sk_b)^2 + C_0(2\gamma + 1)(k - sk_b) + \left[2C_0(k - sk_b) + C_1(k - sk_b)^2 + \frac{2}{3}C_2(k - sk_b)^3\right] \ln|k - sk_b| + i\frac{\pi}{2}s\left[2C_0(k - sk_b) + C_1(k - sk_b)^2 + \frac{2}{3}C_2(k - sk_b)^3\right] \operatorname{sgn}(k - sk_b). \quad (50) $$

Rearranging and differentiating, we obtain, using

$$ \eta_{-1} = 4c_2k_b + 2ic_1k_b^2 - 2c_0k_b^3/3 \quad (51) $$

$$ \eta_0 = 4c_2 + 4ic_1k_b - 2c_0k_b^2 \quad (52) $$

$$ \eta_1 = 2ic_1 - 2c_0k_b \quad (53) $$

$$ \eta_2 = -2c_0/3 \quad (54) $$

$$ \begin{aligned} -D_{k_i} D_{k_j} (k^{-1} \frac{d\hat{v}_1}{dk}) = & \frac{R_{i,j}}{2} + i \frac{\pi}{2} \left[ -4 \frac{c_2}{k_b} \delta(|k| - k_b) \frac{k_i k_j}{k^2} + 2\eta_2 \delta_{i,j} \right. \\ & \left. \qquad \left( \frac{\eta_1}{|k|} - \frac{\eta_{-1}}{|k|^3} \right) \delta_{i,j} + \left( -\frac{\eta_1}{|k|} + 3 \frac{\eta_{-1}}{|k|^3} \right) \frac{k_i k_j}{k^2} \right] \quad (\text{if } |k| \le k_b) \\ & \quad \left( \frac{\eta_1}{|k|} - \frac{\eta_{-1}}{|k|^3} \right) \delta_{i,j} + \left( -\frac{\eta_1}{|k|} + 3 \frac{\eta_{-1}}{|k|^3} \right) \frac{k_i k_j}{k^2} \quad (\text{if } |k| \ge k_b) \end{aligned} \quad (55) $$

with

$$ R_{i,j} = \left[ \frac{2}{3} c_0 \ln \left| \frac{k^2 - k_b^2}{k^2} \right| - \frac{\eta_{-1}}{k_b k^2} + \frac{1}{2} \left( \eta_1 k - \frac{\eta_{-1}}{k} \right) \frac{1}{k^2} \ln \left| \frac{k - k_b}{k + k_b} \right| \right] \delta_{i,j} + \left[ -\frac{1}{2} \left( \eta_1 k - 3 \frac{\eta_{-1}}{k} \right) \frac{1}{k^2} \ln \left| \frac{k - k_b}{k + k_b} \right| + \left( -4c_2 \frac{1}{k^2 - k_b^2} + \frac{3\eta_{-1}}{k_b k^2} \right) \frac{k_i k_j}{k^2} \right] \quad (56) $$
Using Equation (46), we obtain, for the $(i, j)$ component:

$$ 
\begin{aligned}
\mathcal{F}(G_{1Q}(l)\mathbb{Q}_{(i,j)}) &= \left(-\frac{1}{2}\right) \times \left(-D_{k_i} D_{k_j} (k^{-1} \frac{d\hat{v}_1}{dk})\right) \\
&= \left[-\frac{c_0}{6} \ln\left|\frac{k^2 - k_b^2}{k^2}\right| + \frac{\eta_{-1}}{4k_b k^2} - \frac{1}{8}\left(\eta_1 k - \frac{\eta_{-1}}{k}\right) \frac{1}{k^2} \ln\left|\frac{k - k_b}{k + k_b}\right|\right] \delta_{i,j} + \\
&\qquad \left[\frac{1}{8}\left(\eta_1 k - 3\frac{\eta_{-1}}{k}\right) \frac{1}{k^2} \ln\left|\frac{k - k_b}{k + k_b}\right| + \left(c_2 \frac{1}{k^2 - k_b^2} - \frac{3\eta_{-1}}{4k_b k^2}\right) \frac{k_i k_j}{k^2} + \right. \\
&\qquad \left. i \frac{\pi}{2} \left[\frac{c_2}{k_b} \delta(|k| - k_b) \frac{k_i k_j}{k^2} + \left\{ -\frac{\eta_{-1}}{4} \delta_{i,j} \frac{\eta_1}{|k|} \frac{\eta_{-1}}{|k|^3} \right. \right. \\
&\qquad \left. \left. \text{ if } |k| < k_b \right. \right. \\
&\qquad \left. \left. + \frac{\eta_1}{4} \delta_{i,j} + \frac{1}{4} \left(\frac{\eta_1}{|k|} - 3\frac{\eta_{-1}}{|k|^3}\right) \frac{k_i k_j}{k^2} \text{ if } |k| > k_b \right\} \right].
\end{aligned}
\quad (57) $$

From this, using Equations (40)–(43), we can deduce:

$$ 
\begin{aligned}
\mathcal{F}(G_{1I}(l)\mathbb{I} + G_{1Q}(l)\mathbb{Q}) &= \left(a_0 + \frac{c_0}{3}\right) P(k)\mathbb{I} + k_b^{-2} A(k) (a_2\mathbb{I} + c_2\hat{\mathbb{Q}}) \\
&\quad - \frac{ia_1 - a_0k_b}{2} \frac{M(k)}{k}\mathbb{I} + \frac{\eta_{-1}}{8} \left(\frac{2}{k_b} + \frac{M(k)}{k}\right) k^{-2} [\mathbb{I} - 3\hat{\mathbb{Q}}] \\
&\quad + \frac{\eta_1}{8} \frac{M(k)}{k} [\hat{\mathbb{Q}} - \mathbb{I}] .
\end{aligned}
\quad (58) $$

### 4.3.2. Expression in the General Notation

To express Equation (58) in the general notation for the Fourier Transform, let us first define

$$ \tilde{a}_i = \frac{a_i}{k_b^i} $$

$$ \tilde{c}_i = \frac{c_i}{k_b^i} $$

The expression in general notation is obtained from Equation (58) replacing, $a_i$ by $\tilde{a}_i$, $c_i$ by $\tilde{c}_i$, $k_b$ by $k_b/|b|$, and multiplying the whole expression by $c_{a,b}^3$.

## 5. Discussion

### 5.1. Integrability Condition

In order to have a Fourier Transform, the singular part should necessarily be integrable on the unit sphere $S^{n-1}$, which imposes a condition on the zeroth order coefficients only. More precisely, one requires that [24]

$$ \int_{S^{n-1}} \mathbb{G}_s \, dS = 0 . $$

This condition is always satisfied for the odd tensor terms in $\mathbb{G}_s$, but for even tensors, up to order 2, one obtains:

$$ \int_{S^{n-1}} (a_0 \mathbb{I} + c_0 \mathbb{Q}) \, dS = 0 . $$

which imposes

$$ c_0 = -3a_0 , $$

which is satisfied by the Maxwell kernel Equation (4).

### 5.2. At Infinity: Sommerfeld's Rule

The function $r \mapsto E$ should verify some limit condition, which corresponds to scattering by an object in a homogeneous space. Far enough from the scatterer [32], the field $E$ should have
only transverse components, e.g., $Q E(l) \underset{l \to +\infty}{\longrightarrow} 0$. As this should be fulfilled independently of the object shape $\chi$, the kernel should verify $Q G \underset{l \to +\infty}{\longrightarrow} 0$. Noting that $Q$ is a projector, and assuming that the integrability condition, Section 5.1, is satisfied, this imposes

$$a_2 = -c_2$$

Note that this condition is satisfied by the Maxwell kernel Equation (4).

### 5.3. Electromagnetic Scattering

Considering that, for the Maxwell scattering kernel,

$$ 
\begin{aligned}
a_0 &= -1 & c_0 &= 3 \\
a_1 &= ik_b & c_1 &= -3ik_b \\
a_2 &= k_b^2 & c_2 &= -k_b^2
\end{aligned}
$$

one obtains $\eta_{-1} = \eta_1 = 0$.

This leads to the following expression for the FT:

$$ \mathcal{F}(\mathbb{G}_1) = \left( \text{v.p.} \frac{k_b^2}{k^2 - k_b^2} + i \frac{\pi}{2} k_b \delta_{S(0,k_b)}(k) \right) \left[ \mathbb{I} - \hat{\mathbb{Q}} \right], \quad (59) $$

Using the general Fourier notation defined Equation (A1), we obtain:

$$ \mathcal{F}(\mathbb{G}_1) = c_{a,b}^3 A_{a,b}(k) \left[ \mathbb{I} - \hat{\mathbb{Q}} \right], \quad (60) $$

$$ A_{a,b}(k) = \left( \text{v.p.} \frac{k_b^2}{b^2 k^2 - k_b^2} + i \frac{\pi}{2} \frac{k_b}{|b|} \delta_{S(0, k_b / |b|)}(k) \right), \quad (61) $$

where $c_{a,b}^3$ is defined by Equation (16).

Finally, we obtain the Fourier Transform of the LS Equation for scattering, using the general Fourier Transform notations as:

$$ \widehat{\mathbf{E}}_{\text{inc}}(\mathbf{k}, \omega) = \widehat{\mathbf{E}}_{(\mathbf{k}, \omega)} + \left( \frac{|b|}{(2\pi)^{a+1}} \right)^{\frac{3}{2}} \left[ A_{a,b}(k) \left( \hat{\mathbb{Q}} - \mathbb{I} \right) + \hat{\mathbb{Q}} \right] \left( \hat{\chi} \star \hat{\mathbf{E}} \right), \quad (62) $$

with $k = \sqrt{k_x^2 + k_y^2 + k_z^2}$, and $\mathcal{F}(\chi)$ is, e.g., in the case of an isotropic dielectric, the (3D) Fourier Transform of $(x, y, z; \omega) \mapsto \varepsilon(x, y, z; \omega) / \varepsilon_b(\omega) - 1$, where $\varepsilon_b(\omega)$ is the background permittivity. In the anisotropic case, both $\chi$ and $\hat{\chi}$ are $3 \times 3$ matrix.

For a scattering problem, this equation relates the scattered vectorial quantity to the Fourier transform of the object shape, described analytically by its indicator function.

#### 5.3.1. Solution for a Spherical Scatterer, in the Low Contrast Limit

To give an example, we analytically solved Equations (61) and (62). Let us first consider the limit of low index contrast. Note that this corresponds to typical experimental conditions in biophotonics [37,38].

In this case, $\hat{\chi}(\mathbf{k}) \ll 1$ in $\mathbb{R}^3$, and Equation (62) can be inverted easily to give:

$$ \widehat{\mathbf{E}}_{\text{inc}} = \widehat{\mathbf{E}} + \left( \frac{|b|}{(2\pi)^{a+1}} \right)^{\frac{3}{2}} \left[ A_{a,b}(k) \left( \hat{\mathbb{Q}} - \mathbb{I} \right) + \hat{\mathbb{Q}} \right] \left( \hat{\chi} \star \widehat{\mathbf{E}}_{\text{inc}} \right). \quad (63) $$
Considering a plane wave excitation, and using, in the following, the notation $\hat{u}$ for the FT obtained with $a = 1, b = -1$ (see Appendix A),

$$ \hat{\mathbf{E}}_{\text{inc}} = (2\pi)^3 \mathbf{E}_o \delta(\mathbf{k} - \mathbf{k}_b). $$

Then, $\hat{\chi} \star \hat{\mathbf{E}}_{\text{inc}} = (2\pi)^3 \mathbf{E}_o \hat{\chi}(\mathbf{k} - \mathbf{k}_b)$, and the scattered field, defined as $\mathbf{E}_s \equiv \mathbf{E} - \mathbf{E}_{\text{inc}}$ is given by:

$$ \hat{\mathbf{E}}_{s}(\mathbf{k}, \omega) = - \left[ \frac{k_b^2 \hat{\chi}(\mathbf{k} - \mathbf{k}_b)}{k^2 - k_b^2} + i \frac{\pi}{2} k_b \delta(k - k_b) \hat{\chi}(\mathbf{k} - \mathbf{k}_b) \right] (\hat{\mathbb{Q}} - \mathbb{I}) \mathbf{E}_o - \hat{\chi}(\mathbf{k} - \mathbf{k}_b) \hat{\mathbb{Q}} \mathbf{E}_o \quad (64) $$

This equation gives the analytic solution, in Fourier space, for the (vectorial) field scattered by an arbitrarily shaped scatterer in the regime of small index contrast.

Thanks to, e.g., numerical inverse Fourier Transform, one can obtain the 3D field plots in real space, see Figures 1 and 2, using, for a dielectric sphere of radius $\rho$,

$$ \hat{\chi}(\mathbf{k}) \equiv \hat{\chi}(\|\mathbf{k}\|) = \frac{\varepsilon - \varepsilon_b}{\varepsilon_b} j_1(\|\mathbf{k}\|\rho) / (\|\mathbf{k}\|\rho), $$

where $j_1$ is the spherical Bessel function of order 1. The calculation obtained from Equation (64) is in good agreement with Mie theory—see Figures 1 and 2.

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Maps of the real part of the scattered field, from Mie theory—panel (a,c)—and inverse FT of Equation (64)—panel (b,d). A plane wave polarized linearly along $\mathbf{x}$, and propagating along $\mathbf{z}$, towards $z > 0$, is impinging on a sphere of diameter $d = 1.5\lambda$, where $\lambda$ is the vacuum wavelength. The sphere has an index $n_{\text{int}} = 1.35$, and for the surrounding medium $n_{\text{ext}} = 1.33$. Panel (a,b) displays the $\mathbf{x}$-component, panel (c,d) displays the $\mathbf{z}$-component. The black dotted circle figures the sphere boundary.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Real (black) and imaginary (magenta) part of the x component of the scattered field, along the direction of propagation. Parameters are those of Figure 1. Symbols show the results computed with Mie theory, and lines are for the inverse FT of Equation (64).</figcaption></figure>

One notice that the injected field has no z component, and the latter appears in the scattering process, see Figure 1c,d. This effect, though tenuous, is well described by our tensorial approach.

One readily notices that the scattered field is a linear functional of $\chi$, in Equation (64), which is a consequence of the small index contrast hypothesis. In the case of multiple, non-overlapping scatterers, the field would merely be the sum of all individual contributions, stemming from each scatterer.

### 5.3.2. Scattering by an Ensemble of Spheres or Ellipsoids: A Born Expansion

In the case of a large contrast of permittivity between the scatterer and its environment, Equation (62) can be solved using a Born expansion [39]. We proceed here in the Fourier space, rewriting Equation (62) as:

$$ \hat{\mathbf{E}} = \hat{\mathbf{E}}^{(0)} + \hat{\mathbf{E}}^{(1)} + \dots, $$

with

$$ \hat{\mathbf{E}}^{(0)} = \hat{\mathbf{E}}_{\text{inc}}, \text{ and} $$

$$ \hat{\mathbf{E}}^{(j+1)} = -\frac{1}{(2\pi)^3} \left[ A_{1,-1}(k) \left( \hat{\mathbf{Q}} - \mathbf{I} \right) + \hat{\mathbf{Q}} \right] \left( \hat{\chi} \star \hat{\mathbf{E}}^{(j)} \right). \quad (65) $$

As an example, we compute the multiple scattering by an ensemble of identical spheres, in the long wavelength limit. To proceed, we now define

$$ \hat{\chi}(\mathbf{k}) = \frac{\varepsilon - \varepsilon_b}{\varepsilon_b} j_1(\|\mathbf{k}\|\rho) / (\|\mathbf{k}\|\rho) \times \sum_{j=0}^{N} \exp(i\mathbf{k} \cdot \mathbf{R}_j), $$

where the $\mathbf{R}_j$ are vectors of $\mathbb{R}^3$. $\chi$ now describes N non-overlapping spheres, placed at $\mathbf{R}_j$, and of radius $\rho$.
We remark that when several scatterers, described by $(\chi_p)_{p \in [1..P]}$, are excited together (e.g., by a plane wave), one can write $\chi = \sum_p \chi_p$, and express analytically, thanks to Equation (65), all the contributions of the multiple scattering at an order $j$th of the Born expansion.

In the long wavelength limit, where $A_{1,-1}(k) = 0$, the successive Born expansions can be computed iteratively using Fast Fourier Transform. We have used the following algorithm, based on the fact that $\hat{\chi} \star \mathbf{E}^{(j)}$ is the FT of the internal field—the field inside the scatterers—computed at the $j$th order:

(1) Knowing $\mathbf{E}^{(j)}$ in real space, compute the internal field $\chi(\mathbf{r}) \cdot \mathbf{E}^{(j)}(\mathbf{r})$,

(2) Make a 3D FFT to obtain $\hat{\chi} \star \hat{\mathbf{E}}^{(j)}$ in Fourier space,

(3) Multiply the result by $Q$, one obtains $\hat{\mathbf{E}}^{(j+1)}$,

(4) Invert the 3D FFT to obtain the field $\mathbf{E}^{(j+1)}$ in real space, inside *and* outside the scatterers.

(5) Iterate to step (1), with $j \leftarrow j + 1$

One then adds up the orders of Born expansion to obtain the scattered field.

We remark that the first step relies on the internal field only, as $\chi(\mathbf{r}) = 0$ outside the scatterers. Therefore, no specific boundary condition—e.g., perfectly matched layers—is necessary in the real space for this algorithm. However, to avoid the Gibbs effect [40,41], in the implementation of the algorithm, one has to pay attention to the expansion in the Fourier space of the product $\chi(\mathbf{r}) \cdot \mathbf{E}^{(j)}(\mathbf{r})$. As known from previous studies in 1D [2,41] and 2D [42], the unavoidable truncature (due to a finite extent of the numerical Fourier space) can dramatically affect the value of $\hat{\chi} \star \hat{\mathbf{E}}^{(j)}$ that is computed. In the present work, to reconstruct an accurate field profile, we used a specific truncation in the Fourier domain, choosing an apodization function ([43], See Eq. (2)) that permits us to avoid the non-physical Gibbs oscillations. This truncation has been enforced on the Fourier components of $\hat{\chi}$, at step (2) of the algorithm.

Eventually, our results, see Figure 3, are in good agreement with the standard T-matrix method [20] to compute multiple scattering by an ensemble of spheres in a homogeneous medium. Note that only three spheres among the five we considered are visible on the Y-Z plane cut (X = 0), and only one is visible on the X-Z plane cut (Y = 0). This could permit computing the optical forces that drive the dynamics of a trapped ensemble of glass spheres in a vacuum [44–46].

<figure><img src="image_7.png" /><figcaption>**Figure 3.** Map of the intensity of the x component of the scattered field, in the (XZ) and (YZ) planes, for a cluster of five spheres of glass ($\epsilon = 2.1$) in air ($\epsilon_b = 1$), placed at [0,0,0]; [0,R,2R]; [0.2R,−2R,R]; [2R,R,−2R]; [−1.5R,−1.3R,2R], excited by a plane wave polarized linearly along Z, propagating along X, of unitary amplitude, with $R = \lambda/100$. Exact result, from [20] (panel (a,c)) and results obtained from the fourth order Born expansion of Equation (62) (panel (b,d)) are displayed.</figcaption></figure>
Finally, using the same method, we are able to compute the field scattered by an ensemble of randomly oriented ellipsoid of semi-axes (a, b, c), whose shape can be described analytically in the Fourier space by:

$$ \hat{\chi}(\mathbf{k}) = \sum_{j=0}^{N} \exp(i\mathbf{k} \cdot \mathbf{R}_j) \frac{\varepsilon - \varepsilon_b}{\varepsilon_b} \frac{j_1(\sqrt{\mathbf{k} \cdot \mathbf{M}(\phi_j, \theta_j, \psi_j) \mathbf{k}^t})}{\sqrt{\mathbf{k} \cdot \mathbf{M}(\phi_j, \theta_j, \psi_j) \mathbf{k}^t}}, \quad (66) $$

where $M(\phi_j, \theta_j, \psi_j)$ is the descriptor for the ellipsoid, rotated by the Euler angles $\phi_j$ around X axis, $\theta_j$ around Y axis and $\psi_j$ around Z axis:

$$ M(\phi_j, \theta_j, \psi_j) = R(\phi_j, \theta_j, \psi_j) \begin{pmatrix} a^2 & 0 & 0 \\ 0 & b^2 & 0 \\ 0 & 0 & c^2 \end{pmatrix} R(\phi_j, \theta_j, \psi_j)^{-1}, $$

$R(\phi_j, \theta_j, \psi_j)$ being the associated 3 × 3 Euler rotation matrix.

The result is presented in Figure 4. One can see that the interaction between neighboring scatterers, depending on their orientation, can deeply modify the field locally, e.g., enhance some components of the scattered field that were not present in the injected field.

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Map of the intensity of two components of the scattered field, in the (YZ) plane, for a cluster of 20 randomly placed and oriented ellipsoid of glass ($\varepsilon = 2.1$) in air ($\varepsilon_b = 1$), whose semi-lengths are $0.9R$, $0.5R$, and $0.3R$, with $R = \lambda/100$. The scatterers are excited by a plane wave polarized linearly along Z, propagating along X, of unitary amplitude. The Z-component (panel (a)) and the Y-component (panel (b)) are displayed.</figcaption></figure>

## 6. Conclusions

In summary, we studied the Fourier Transformation of the Lippmann-Schwinger equation for electromagnetic scattering. Dealing separately with the singular and the regular part of the equation, we obtained an analytic expression where the role of the scatterer shape is clearly visible. This does not require a periodic geometry, nor the use of perfectly matched layers to adapt the boundary conditions.

We have provided a few numerical examples, using the factorization in Fourier space to solve the equation for simple systems, with low index contrasts, or in the long wavelength limit, in a configuration where multiple scattering occurs.

Perspectives include the possibility to solve, in Fourier space, the Lippmann-Schwinger equation associated with scattering problems in more general cases, taking advantage of recent works to compute fast Fourier Transforms [47]. This could permit, e.g., defining the Quasi-Normal Modes [18,48,49] of a dispersive scattering system. Equation (62) could help, e.g., to obtain a relationship between the scatterer shape and its eigenfunctions (the QNM modes), in particular for the static modes [50], for which $k_b = 0$. Field computation could also be carried out in lower dimensional systems, such as structured fibers [51]. Finally, this
approach could help us to understand the multiple scattering in an ensemble of particles trapped by light beams using an analytical formalism.

**Author Contributions:** Conceptualization, F.G.; Methodology, F.G., V.R. and M.P.; Software, M.P.; Validation, F.G. and M.P.; Formal analysis, F.G. and V.R.; Investigation, F.G., V.R. and M.P.; Writing—original draft, F.G. and V.R.; Writing—review & editing, F.G. and M.P.; Supervision, F.G. and M.P.; Project administration, M.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Agence Nationale pour la Recherche (Grant No ANR-21-CE30-0006 and No ANR-22-CE24-0010).

**Data Availability Statement:** Data are available upon reasonable request.

**Acknowledgments:** M.P. acknowledges Suzanne Bessus, Josiane Parzych, Valérie Thouard and Sophie Sauvant for administrative support.

**Conflicts of Interest:** The authors declare no conflict of interest.

# Appendix A. Fourier Transform and Notations

## Appendix A.1. Fourier Transform of a Function

There are a lot of definitions of the Fourier transform. Let us present a convenient notation to encompass all the cases once for all (with $f,g: \mathbb{R}^n \to \mathbb{R}$, $\mathcal{F}$ representing the direct Fourier transform). Let us denote $\langle f,g \rangle$ the canonical scalar product, either on $\mathbb{R}^n$, or (unambiguously) on the function space $\mathbb{R}^n \to \mathbb{R}$. In the following, the fact that the electromagnetic field is a complex valued function will be handled by splitting real and imaginary part, in the calculations:

$$ \mathcal{F}^{a,b} f(k) := \left( \frac{|b|}{(2\pi)^{1-a}} \right)^{n/2} \int_{\mathbb{R}^n} f(t) e^{ib\langle k, t \rangle} \, dt, \quad (A1) $$

$$ (\mathcal{F}^{a,b})^{-1} g(t) := \left( \frac{|b|}{(2\pi)^{1+a}} \right)^{n/2} \int_{\mathbb{R}^n} g(k) e^{-ib\langle t, k \rangle} \, dk, \quad (A2) $$

These simple results are gathered here withoug proof. For more details, see [31].

To jump from a convention to another, we have the simple following corresponding formula:

$$ \mathcal{F}^{a,b} f(k) = \left( \frac{|b/b'|}{(2\pi)^{a'-a}} \right)^{n/2} \mathcal{F}^{a',b'} f\left(\frac{b}{b'}k\right). \quad (A3) $$

## Appendix A.2. Fourier Transform of a Tempered Distribution

Let us recall that if $T \in S'(\mathbb{R}^n)$ is a tempered distribution, we define the Fourier transform $\mathcal{F}^{a,b} T$ by

$$ \langle \mathcal{F}^{a,b} T, \varphi \rangle := \langle T, \mathcal{F}^{a,b} \varphi \rangle, \quad (A4) $$

where,

$$ \langle f,g \rangle := \int_{\mathbb{R}^n} f(t)g(t) \, dt. \quad (A5) $$

One can generalize the “jump” formula Equation (A3) to distribution, by:

$$ \mathcal{F}^{a,b} T = \left( \frac{|b/b'|}{(2\pi)^{a'-a}} \right)^{n/2} m_{\frac{b}{b'}} \mathcal{F}^{a',b'} T, \quad (A6) $$

where, for all $\alpha \neq 0$,

$$ \langle m_{\alpha} T, \varphi \rangle := \frac{1}{|\alpha|^n} \langle T, m_{\frac{1}{\alpha}} \varphi \rangle \quad (A7) $$
### Appendix A.3. Fourier Transform of a Convolution

The FT of convolution can be expressed as:

$$ \mathcal{F}^{a,b}(f \star g) = \left( \frac{(2\pi)^{1-a}}{|b|} \right)^{\frac{n}{2}} \left( \mathcal{F}^{a,b} f \right) \cdot \left( \mathcal{F}^{a,b} g \right), \qquad (A8) $$

### Appendix A.4. Fourier Transform of a Product

The well known formula that links the FT of a product and the convolution of two FTs can be casted under the general notation, see Equations (A1) and (A2). The result reads:

$$ \mathcal{F}^{a,b}(f \cdot g) = \left( \frac{|b|}{(2\pi)^{a+1}} \right)^{\frac{n}{2}} \left( \mathcal{F}^{a,b} f \right) \star \left( \mathcal{F}^{a,b} g \right), \qquad (A9) $$

### Appendix A.5. Alternative Notation Used in This Article

These general notations, in spite of being really useful to bring altogether the different Fourier transform conventions, can be a bit heavy while conducting computations. We will often use, for a function (or a distribution) $u$ the simplified notation $\mathcal{F}(u)$, or even $\hat{u}$, to symbolize $\mathcal{F}^{1,-1}(u)$.

## Appendix B. Lemma on Spherical Properties

### Appendix B.1. Spherical Properties

Here, we give some detail on the notation we use in Appendix B.2.

In the multidimensional case, we have

$$ \int_{\mathbb{S}^{n-1}} f(\theta) \, d\sigma(\theta) = \int_{\varphi_1=0}^{\pi} \dots \int_{\varphi_{n-2}=0}^{\pi} \int_{\varphi_{n-1}=0}^{2\pi} f(\theta(\varphi)) J(n, \varphi) \, d\varphi_{n-1} \dots d\varphi_1 \qquad (A10) $$

where

$$ \begin{aligned} \theta_1 &= \cos(\varphi_1), \\ \theta_2 &= \sin(\varphi_1) \cos(\varphi_2), \\ \theta_3 &= \sin(\varphi_1) \sin(\varphi_2) \cos(\varphi_3), \\ \vdots \\ \theta_{n-1} &= \sin(\varphi_1) \dots \sin(\varphi_{n-2}) \cos(\varphi_{n-1}), \\ \theta_n &= \sin(\varphi_1) \dots \sin(\varphi_{n-2}) \sin(\varphi_{n-1}), \end{aligned} $$

and

$$ J(n, \varphi) = (\sin(\varphi_1))^{n-2} \dots (\sin(\varphi_{n-3}))^2 \sin(\varphi_{n-2}). \qquad (A11) $$

Note that $\varphi_1, \varphi_2, \dots, \varphi_{n-2} \in [0, \pi]$ and $\varphi_{n-1} \in [0, 2\pi]$.

### Appendix B.2. Detailed Computation of the Integrals for $\mathcal{G}_s$

First, we can manage to obtain $\int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_2^2 d\theta$ from $\int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_1^2 d\theta$ and $\int_{\mathbb{S}^{n-1}} \log |\theta_1| d\theta$.

Using the notation defined in Appendix B.1, one obtains:

$$ \begin{aligned} \int g(\varphi_1) \cos^2(\varphi_2) J(n, \varphi) \, d\varphi &= \int g(\varphi_1) J(n, \varphi) \, d\varphi - \int g(\varphi_1) \sin^2(\varphi_2) J(n, \varphi) \, d\varphi \\ &= \int g(\varphi_1) J(n, \varphi) \, d\varphi - \frac{n-2}{n-1} \int g(\varphi_1) J(n, \varphi) \, d\varphi \\ &= \frac{1}{n-1} \int g(\varphi_1) J(n, \varphi) \, d\varphi \end{aligned} $$
Using Fubini and the properties of Wallis integrals: $W_n = \int_0^{\frac{\pi}{2}} \sin^n(x) dx$ that obey the well-known relation $(n-1)W_{n-2} = nW_n$, for all $n \ge 2$. (here $n \ge 4$),

$$ 
\begin{aligned}
\int_0^{\pi} \sin^2(\varphi_2) \sin(\varphi_2)^{n-3} d\varphi_2 &= 2W_{n-1} \\
&= \frac{n-2}{n-1} 2W_{n-3} \\
&= \frac{n-2}{n-1} \int_0^{\pi} \sin(\varphi_2)^{n-3} d\varphi_2
\end{aligned}
$$

So

$$ 
\begin{aligned}
\int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_2^2 d\theta &= \int \log |\cos(\varphi_1)| \sin^2(\varphi_1) \cos^2(\varphi_2) J(n, \varphi) d\varphi \\
&= \frac{1}{n-1} \int \log |\cos(\varphi_1)| \sin^2(\varphi_1) J(n, \varphi) d\varphi \\
&= \frac{1}{n-1} \left( \int_{\mathbb{S}^{n-1}} \log |\theta_1| d\theta - \int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_1^2 d\theta \right)
\end{aligned}
$$

To compute these integrals, we use ([34], p. 442), that links the integral of $\theta \mapsto K(x \cdot \theta)$ on the $\mathbb{S}^{n-1}$ sphere to the following integral on the $[0 \ 1]$ interval:

$$ 
\begin{aligned}
\int_0^1 s^\alpha (1-s^2)^\beta \log(s) ds &= \frac{\partial}{\partial \alpha} \int_0^1 s^\alpha (1-s^2)^\beta ds \\
&= \frac{\partial}{2\partial \alpha} \int_0^1 u^{(\alpha-1)/2} (1-u)^\beta du \\
&= \frac{\partial}{2\partial \alpha} B\left(\frac{\alpha+1}{2}, \beta+1\right) \\
&= \frac{1}{4} B\left(\frac{\alpha+1}{2}, \beta+1\right) \left(\psi\left(\frac{\alpha+1}{2}\right) - \psi\left(\frac{\alpha+1}{2} + \beta + 1\right)\right)
\end{aligned}
$$

where $\psi(x) := \frac{\Gamma'(x)}{\Gamma(x)}$ is the digamma function, and $B$ is defined by $B(x, y) = \frac{\Gamma(x)\Gamma(y)}{\Gamma(x+y)}$, and verifies $\frac{\partial}{\partial x}B(x, y) = B(x, y)(\psi(x) - \psi(x+y))$. Then, using the relations $\Gamma(z+1) = z\Gamma(z)$, $\Gamma(\frac{1}{2}) = \sqrt{\pi}$, $\psi(z+1) = \psi(z) + \frac{1}{z}$, $\psi(\frac{1}{2}) = -2\log(2) - \gamma$, we obtain:

$$ \int_{\mathbb{S}^{n-1}} \log |\theta_1| d\theta = - \frac{\pi^{n/2}}{\Gamma(n/2)} \left( \psi(n/2) + \gamma + 2\log(2) \right) \quad (A12) $$

and

$$ \int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_1^2 d\theta = - \frac{\pi^{n/2}}{n\Gamma(n/2)} \left( \psi(n/2) + \frac{2}{n} + \gamma + 2\log(2) - 2 \right). \quad (A13) $$

We also obtain

$$ 
\begin{aligned}
\int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_2^2 d\theta &= \frac{1}{n-1} \left( \int_{\mathbb{S}^{n-1}} \log |\theta_1| d\theta - \int_{\mathbb{S}^{n-1}} \log |\theta_1| \theta_1^2 d\theta \right) \\
&= - \frac{\pi^{n/2}}{n\Gamma(n/2)} \left( \psi(n/2) + \frac{2}{n} + \gamma + 2\log(2) \right).
\end{aligned}
\quad (A14) $$

Using the same method, one computes without difficulty the integrals for higher order tensor—the results are gathered Table 2.

## References

1. Goodman, J. *Introduction to Fourier Optics*, 2nd ed.; Mc Graw-Hill: New York, NY, USA, 1996.

2. Lalanne, P.; Morris, G.M. Highly improved convergence of the coupled-wave method for TM polarization. *J. Opt. Soc. Am. A* **1996**, *13*, 779. [CrossRef]

3. Granet, G.; Guizal, B. Efficient implementation of the coupled-wave method for metallic lamellar gratings in TM polarization. *J. Opt. Soc. Am. A* **1996**, *13*, 1019–1023. [CrossRef]
4. Knop, K. Rigorous diffraction theory for transmission phase gratings with deep rectangular grooves. *J. Opt. Soc. Am.* **1978**, *68*, 1206–1210. [CrossRef]

5. Popov, E.; Nevière, M. Maxwell equations in Fourier space: Fast-converging formulation for diffraction by arbitrary shaped, periodic, anisotropic media. *J. Opt. Soc. Am. A* **2001**, *18*, 2886–2894. [CrossRef]

6. Bertrand, A.; Dumur, F.; Mruczkiewicz, M.; Perrin, M.; Lartigau-Dagron, C.; Bousquet, A.; Vignau, L.; Billon, L.; Fasquel, S. Bottom-up honeycomb top layer for light outcoupling enhancement in blue organic light emitting diodes. *Org. Electron.* **2018**, *52*, 222–229. [CrossRef]

7. Dumur, F.; Reculusa, S.; Mruczkiewicz, M.; Perrin, M.; Vignau, L.; Fasquel, S. Multilayer Langmuir-Blodgett films as diffractive external 3D photonic crystal in blue OLEDs. *Opt. Express* **2016**, *24*, 27184–27198. [CrossRef]

8. Bachelet, A.; Fasquel, S.; Rampnoux, J.M.; Jonusauskas, G.; Takimiya, K.; Hirsch, L.; Perrin, M.; Abbas, M. Wide Range Color Tuning in Single Emissive Layer Organic Light Emitting Transistors. *ACS Photonics* **2023**, *10*, 2793–2798. [CrossRef]

9. Arnoldus, H.F. Representation of the near-field, middle-field, far-field electromagnetic Green's functions in reciprocal space. *J. Opt. Soc. Am. B* **2001**, *18*, 547–555. [CrossRef]

10. Sheppard, C.J.R.; Lin, J.; Kou, S.S. Rayleigh-Sommerfeld diffraction formula in k space. *J. Opt. Soc. Am. A* **2013**, *30*, 1180–1183. [CrossRef]

11. Sheppard, C.J.R.; Kou, S.S.; Lin, J. The Green-function transform and wave propagation. *Front. Phys.* **2014**, *2*, 67. [CrossRef]

12. Schmalz, J.A.; Schmalz, G.; Gureyev, T.E.; Pavlov, K.M. On the derivation of the Green's function for the Helmholtz equation using generalized functions. *Am. J. Phys.* **2010**, *78*, 181. [CrossRef]

13. Acquista, C. Light scattering by tenous particles: A generalization of the Rayleigh-Gans-Rocard approach. *Appl. Opt.* **1976**, *15*, 2932–2936. [CrossRef] [PubMed]

14. Gallatin, G.M. Fourier, gauss, fraunhofer, porod and the shape from moments problem. *J. Math. Phys.* **2012**, *53*, 013509. [CrossRef]

15. Ciccariello, S.; Schneider, J.-M.; Schonfeld, B.; Kostorz, G. Generalization of Porod's law of small-angle scattering to anisotropic samples. *Euro Phys. Lett.* **2000**, *50*, 601–607. [CrossRef]

16. Vynck, K.; Pierrat, R.; Carminati, R. Polarization and spatial coherence of electromagnetic waves in uncorrelated disordered media. *Phys. Rev. A* **2014**, *89*, 013842. [CrossRef]

17. Kahnert, F.M. Numerical methods in electromagnetic scattering theory. *J. Quant. Spect. Rad. Transf.* **2003**, 79–80, 775–824. [CrossRef]

18. Bai, Q.; Perrin, M.; Sauvan, C.; Hugonin, J.; Lalanne, P. Efficient and intuitive method for the analysis of light scattering by a resonant nanostructure. *Opt. Express* **2013**, *21*, 27371–27382. [CrossRef]

19. Waterman, P.C. Matrix formulation of electromagnetic scattering. *Proc. IEEE* **1965**, *53*, 805–812. [CrossRef]

20. Egel, A.; Pattelli, L.; Mazzamuto, G.; Wiersma, D.S.; Lemmer, U. CELES: CUDA-accelerated simulation of electromagnetic scattering by large ensembles of spheres. *J. Quant. Spect. Rad. Transf.* **2017**, *199*, 103–110. [CrossRef]

21. Mishchenko, M.I.; Travis, L.D.; Mackowski, D.W. T-matrix computations of light scattering by nonspherical particles: A review. *J. Quant. Spectr. Rad. Transf.* **1996**, *55*, 535–575. [CrossRef]

22. Hu, S.; Liu, L.; Gao, T.; Zeng, Q. Design and Validation of the Invariant Imbedded T-Matrix Scattering Model for Atmospheric Particles with Arbitrary Shapes. *Appl. Sci.* **2019**, *9*, 4423. [CrossRef]

23. Schebarchov, D.; Ru, E.C.L.; Grand, J.; Auguié, B. Mind the gap: Testing the Rayleigh hypothesis in T-matrix calculations with adjacent spheroids. *Opt. Express* **2019**, *27*, 35750–35760. [CrossRef]

24. Budko, N.V.; Samokhin, A.B. Spectrum of the volume integral operator of electromagnetic scattering. *SIAM J. Sci. Comput.* **2006**, *28*, 682–700. [CrossRef]

25. Zouros, G.P.; Budko, N.V. Transverse electric scattering on inhomogeneous objects: Spectrum of integral operator and preconditioning. *SIAM J. Sci. Comput.* **2012**, *34*, B226–B246. [CrossRef]

26. Baddour, N. Operational and convolution properties of three-dimensional Fourier transforms in spherical polar coordinates. *J. Opt. Soc. Am. A* *Opt. Image Sci. Vis.* **2010**, *27*, 2144–2155. [CrossRef] [PubMed]

27. Reed, M.; Simon, B. *Methods of Modern Mathematical Physics, Volume 3 of Scattering Theory*; Academic Press: Cambridge, MA, USA, 1979.

28. Samokhin, A.B.; Shestopalov, Y.V.; Shestopalov, I.U.V. *Integral Equations and Iteration Methods in Electromagnetic Scattering*; VSP: Zeist, The Netherlands, 2001.

29. Perrin, M.; Gruy, F. Explicit calculation of singular integrals of tensorial polyadic kernels. *Q. Appl. Math.* **2023**, *81*, 65–86. [CrossRef]

30. Grafakos, L.; Teschl, G. On Fourier transforms of radial functions and distributions. *J. Fourier Anal. Appl.* **2013**, *19*, 167–179. [CrossRef]

31. Gruy, F.; Perrin, M.; Rabiet, V. Fourier Transform of the Lippmann-Schwinger Equation for 3D Vectorial Electromagnetic Scattering: A Direct Relationship between Fields and Shape. Available online: https://hal.archives-ouvertes.fr/hal-03043716 (accessed on 7 December 2020).

32. Schot, S.H. Eighty years of Sommerfeld's radiation condition. *Hist. Math.* **1992**, *19*, 385–401 [CrossRef]

33. Novotny, L.; Hecht, B. *Principles of Nano-Optics*; Cambridge University Press: Cambridge, UK, 2006.

34. Grafakos, L. Classical Fourier Analysis. In *Graduate Texts in Mathematics*, 3rd ed.; Springer: New York, NY, USA, 2014.

35. Rabiet, V. Fourier Transform of a Class of Radial and Semi-Radial Functions. Available online: https://hal.archives-ouvertes.fr/hal-02494837 (accessed on 29 February 2020).
36. Matiur, R. *Applications of Fourier Transforms to Generalized Functions*; WIT Press: Billerica, MA, USA, 2011.

37. Delabre, U.; Feld, K.; Crespo, E.; Whyte, G.; Sykes, C.; Seifert, U.; Guck, J. Deformation of phospholipid vesicles in an optical stretcher. *Soft Matter* **2015**, *11*, 6075–6088. [CrossRef]

38. Crouzil, T.; Perrin, M. Dynamics of a chain of optically coupled micro droplets. *J. Eur. Opt. Soc.-Rapid Publ.* **2013**, *8*, 13079-1. [CrossRef]

39. Kilgore, K.; Moskow, S.; Schotland, J.C. Convergence of the Born and inverse Born series for electromagnetic scattering. *Appl. Anal.* **2017**, *96*, 1737–1748. [CrossRef]

40. Villeneuve, P.R.; Piché, M. Photonic bandgaps: What is the best numerical representation of periodic structures? *J. Mod. Opt.* **1994**, *41*, 241–256. [CrossRef]

41. Li, L. Use of Fourier series in the analysis of discontinuous periodic structures. *J. Opt. Soc. Am. A* **1996**, *13*, 1870–1876. [CrossRef]

42. Popov, E.; Nevière, M.; Bonod, N. Factorization of products of discontinuous functions applied to Fourier–Bessel basis. *J. Opt. Soc. Am. A* **2004**, *21*, 46–52. [CrossRef]

43. Parker, K.J. Apodization and Windowing Functions. *IEEE Trans. Ultrason. Ferro. Freq. Control* **2013**, *60*, 1263–1271. [CrossRef] [PubMed]

44. Amarouchene, Y.; Mangeat, M.; Vidal-Montes, B.; Ondic, L.; Guerin, T.; Dean, D.S.; Louyer, Y. Nonequilibrium Dynamics Induced by Scattering Forces for Optically Trapped Nanoparticles in Strongly Inertial Regimes. *Phys. Rev. Lett.* **2019**, *122*, 183901. [CrossRef]

45. Bellando, L.; Kleine, M.; Amarouchene, Y.; Perrin, M.; Louyer, L. Giant Diffusion of Nanomechanical Rotors in a Tilted Washboard Potential. *Phys. Rev. Lett.* **2022**, *129*, 023602. [CrossRef] [PubMed]

46. Hüpfl, J.; Bachelard, N.; Kaczvinszki, M.; Horodynski, M.; Kühmayer, M.; Rotter, S. Optimal Cooling of Multiple Levitated Particles through Far-Field Wavefront Shaping. *Phys. Rev. Lett.* **2023**, *130*, 083203. [CrossRef]

47. Leitersdorf, O.; Boneh, Y.; Gazit, G.; Ronen, R.; Kvatinsky, S. FourierPIM: High-throughput in-memory Fast Fourier Transform and polynomial multiplication. *Mem.-Mater. Devices Circuits Syst.* **2023**, *4*, 100034. [CrossRef]

48. Perrin, M. Eigen-energy effects and non-orthogonality in the quasi-normal mode expansion of maxwell equations. *Opt. Express* **2016**, *24*, 27137–27151. [CrossRef] [PubMed]

49. Yan, W.; Faggiani, R.; Lalanne, P. Rigorous modal analysis of plasmonic nanoresonators. *Phys. Rev. B* **2018**, *97*, 205422. [CrossRef]

50. Sauvan, C. Quasinormal modes expansions for nanoresonators made of absorbing dielectric materials: Study of the role of static modes. *Opt. Express* **2021**, *29*, 8268–8282. [CrossRef] [PubMed]

51. Bouwmans, G.; Pureur, V.; Betourne, A.; Quiquempois, Y.; Perrin, M.; Bigot, L.; Douay, M. Progress in solid core photonic bandgap fibers. *Opt. Quantum Electron.* **2007**, *39*, 949–961. [CrossRef]

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.