# Eliminating blowing-ups and evanescent waves when using the finite series technique in evaluating beam shape coefficients for some T-matrix approaches, with the example of Gaussian beams

Gérard Gouesbet, Jianqi Shen, Leonardo André Ambrosio

► To cite this version:

Gérard Gouesbet, Jianqi Shen, Leonardo André Ambrosio. Eliminating blowing-ups and evanescent waves when using the finite series technique in evaluating beam shape coefficients for some T-matrix approaches, with the example of Gaussian beams. Journal of Quantitative Spectroscopy and Radiative Transfer, 2024, 330, pp.109212. 10.1016/j.jqsrt.2024.109212. hal-04829046

HAL Id: hal-04829046

https://hal.science/hal-04829046v1

Submitted on 23 Jan 2026

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 CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

# Eliminating blowing-ups and evanescent waves when using the finite series technique in evaluating beam shape coefficients for some T-matrix approaches, with the example of Gaussian beams

Gérard Gouesbet<sup>a,*</sup>, Jianqi Shen<sup>b</sup>, Leonardo André Ambrosio<sup>c</sup>

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

## ARTICLE INFO

**Keywords:**

Generalized Lorenz-Mie theory

T-matrix

Beam shape coefficients

Finite series

Blowing-ups

Evanescent waves

## ABSTRACT

When evaluating beam shape coefficients which encode the description of laser beams, for use in some T-matrix approaches such as generalized Lorenz-Mie theory or Extended Boundary Condition Method for structured beams, by using finite series, blowing-ups are observed. When numerical inaccuracies are ruled out, it has been firmly demonstrated that such blowing-ups correspond to genuine physical phenomena, namely they describe evanescent waves. We propose a method to eliminate these blowing-ups and the corresponding evanescent waves (at least most of them).

## 1. Introduction

In some T-matrix approaches with arbitrary shaped beam illumination, specifically in Generalized Lorenz-Mie Theory (GLMT), e.g. [1–3], or in Extended Boundary Condition Method (EBCM) [4–6] when applied to structured beams [7,8], illuminating beams may be expressed by an expansion over Vector Spherical Wave Functions (VSWFs) and the expansion coefficients are then expressed in terms of Beam Shape Coefficients (BSCs), e.g. [9–11] and references therein for reviews of these theories and applications.

The BSCs are in general constituted from two sets, and are denoted most usually as $g_{n,TM}^m$ and $g_{n,T_E}^m$ with $n$ ranging from 1 to $\infty$, $m$ from $(-n)$ to $(+n)$ and in which *TM* stands for “Transverse Magnetic” while *TE* stands for “Transverse Electric”. There exist several methods to evaluate the BSCs, including (i) the original quadrature method, e.g. [12,13], (ii) the localized approximations which have a long story with a review in [14] to be completed by [15–17], and by a series of paper devoted to the limitations of localized approximations when dealing with beams exhibiting an axicon angle and/or a topological charge, e.g. [18–25], (iii) the use of an Angular Spectrum Decomposition (ASD), e.g. [26] for a review and [27–38], for recent implementations, and (iv) the recently introduced *R*-quadrature method [38–42].

In the present paper, we are concerned with another technique, named the finite series technique. This technique has been introduced several decades ago [43,44] and has been essentially forgotten for a

long time, particularly due to the efficiency of localized approximations in terms of formal flexibility and computing speed. It has recently been set again on the stage due to the above mentioned limitations of localized approximations for beams exhibiting an axicon angle and/or a topological charge, e.g. [45–53].

It has soon been observed that the use of finite series generates blowing-ups in the values of the high-order BSCs, e.g. section 6.7 in [1]. These blowing-ups were attributed to a genuine phenomenon but without any physical explanation of their origins. The renewal of the finite series technique set again the blowing-up issue on the stage. It has then been firmly demonstrated by analytical means that indeed, once numerical inaccuracies are ruled out, blowing-ups are genuine phenomena and, furthermore, it has been demonstrated that they are a signature of the existence of evanescent waves [54,55], and see particularly Eqs. (48)–(49) in [56]. Furthermore, in the case of Gaussian beams characterized by a beam confinement factor *s*, various criteria to evaluate the critical value of the order of BSCs have been established [54]. These criteria will be discussed as well in the present paper.

Specifically, the present paper presents a method to get rid of these blowing-ups and of the associated evanescent waves (at least most of them). It is interesting to note that the method relies on a mathematical step which, in utmost rigor, is not allowed but which nevertheless provides a correct solution to the issue at stake.

<footnote><sup>a</sup> CORIA-UMR 6614- Normandie Université, CNRS-Université et INSA de Rouen, Campus Universitaire du Madrillet, 76800, Saint-Etienne-du Rouvray, France</footnote>

<footnote><sup>b</sup> College of Science, University of Shanghai for Science and Technology, 516, Jungong Road, Shanghai, 200093, China</footnote>

<footnote><sup>c</sup> Department of Electrical and Computer Engineering, São Carlos School of Engineering, University of São Paulo, 400 Trabalhador São-carlene Ave., São Paulo, SP 13566-590, Brazil</footnote>

<footnote>* Corresponding author.
E-mail address: gouesbet@coria.fr (G. Gouesbet).</footnote>
The paper is organized as follows. Section 2 provides a background and poses a problem to be solved in which BSCs have to be evaluated. Section 3 deals with a first method to the evaluation of these BSCs. Section 4 deals with a second method to the evaluation of the same BSCs. In section 5, we provide a comparison of the results obtained in sections 3 and 4, exhibit the fact that they disagree, one exhibiting the evanescent wave behavior and the other getting rid of them, and explain the reason of this difference of behaviors. Section 6 is a conclusion.

## 2. Background to a problem to be solved

We recall that Gaussian beams may be described by using a scheme of successive approximations due to Davis [57–59]. In this scheme, using Cartesian coordinates (x, y, z), a x-polarized potential vector is introduced for a beam propagating along the z-direction. The nonzero x-component of the potential vector is expressed in terms of a potential function $\psi$ which is expressed by an infinite series reading as:

$$ \psi = \sum_{j=0}^{\infty} s^{2j} \psi_{2j} = \psi_0 + s^2 \psi_2 + s^4 \psi_4 + \dots \quad (1) $$

in which s is a small parameter, called the beam confinement factor (or beam shape factor) defining the focusing of the beam, reading as $1/(k w_0)$ in which k is the wavenumber and $w_0$ the beam waist radius of the beam. Eq. (1) allows one to introduce successive approximations called the first-order, third-order, ... approximations. Once an approximation is obtained, the corresponding electric and magnetic fields, E and H respectively, may be obtained, for instance by using the Lorenz gauge.

The first-order approximation is called the L-approximation ("L" standing for lowest) in which the nonzero field components are $E_x$, $E_z$, $H_y$ and $H_z$ given by Eqs. (4.24)–(4.27) in [1]. Neglecting the components $E_z$ and $H_z$, we obtain a still simpler approximation called the L⁻ approximation in [1], section 4.1.6. In this approximation which is to be used in the present paper, the electromagnetic fields read as [1]:

$$ E(\mathbf{r}) = e_x E_0 Q \exp(-Q \rho^2 / w_0^2 - i k z) \quad (2) $$

$$ H(\mathbf{r}) = e_y H_0 Q \exp(-Q \rho^2 / w_0^2 - i k z) \quad (3) $$

in which $e_x$ and $e_y$ are unit vectors in the x- and y-directions respectively, $E_0$ and $H_0$ are field strengths, $\rho^2 = x^2 + y^2$, and $Q$, which is defined in a way slightly differently from in [1], reads as:

$$ Q = \frac{1}{1 - 2iz/(k w_0^2)} \quad (4) $$

The radial components of the electric and magnetic fields may be obtained in spherical coordinates (r, θ, φ) from Eqs. (3.10) and (3.19) in [1]. Together with Eqs. (2) and (3), we then obtain:

$$ \begin{aligned} E_r(\mathbf{r}) &= E_0 Q \exp(-Q r^2 \sin^2 \theta / w_0^2 - i k r \cos \theta) \sin \theta \cos \varphi \\ &= \frac{E_0}{r} \sum_{n=1}^{\infty} \sum_{m=-n}^{+n} c_n^{\text{pw}} n(n+1) g_{n,TM}^m j_n(kr) P_n^{|m|}(\cos \theta) \exp(im\varphi) \end{aligned} \quad (5) $$

$$ \begin{aligned} H_r(\mathbf{r}) &= H_0 Q \exp(-Q r^2 \sin^2 \theta / w_0^2 - i k r \cos \theta) \sin \theta \sin \varphi \\ &= \frac{H_0}{r} \sum_{n=1}^{\infty} \sum_{m=-n}^{+n} c_n^{\text{pw}} n(n+1) g_{n,TE}^m j_n(kr) P_n^{|m|}(\cos \theta) \exp(im\varphi) \end{aligned} \quad (6) $$

in which $j_n(kr)$ are spherical Bessel functions of the first kind, $P_n^{|m|}(\cos \theta)$ are associated Legendre functions defined according to Hobson's notation [60], and $c_n^{\text{pw}}$, with "pw" standing for "plane wave", are plane wave coefficients given by Eq. (3.3) in [1], reading as:

$$ c_n^{\text{pw}} = \frac{1}{k} (-i)^{n+1} \frac{2n+1}{n(n+1)} \quad (7) $$

Introducing $R = kr$ and rearranging, Eqs. (5)–(6) may be rewritten as:

$$ \begin{aligned} & RQ \exp(-Qs^2 R^2 \sin^2 \theta - iR \cos \theta) \sin \theta \cos \varphi \\ &= \sum_{n=1}^{\infty} \sum_{m=-n}^{+n} (-i)^{n+1} (2n+1) g_{n,TM}^m j_n(R) P_n^{|m|}(\cos \theta) \exp(im\varphi) \end{aligned} \quad (8) $$

$$ \begin{aligned} & RQ \exp(-Qs^2 R^2 \sin^2 \theta - iR \cos \theta) \sin \theta \sin \varphi \\ &= \sum_{n=1}^{\infty} \sum_{m=-n}^{+n} (-i)^{n+1} (2n+1) g_{n,TE}^m j_n(R) P_n^{|m|}(\cos \theta) \exp(im\varphi) \end{aligned} \quad (9) $$

We then multiply both sides by $\exp(-im'\varphi)$ and carry out the integration over the interval $[0, 2\pi]$, using classical expressions reading as:

$$ \int_0^{2\pi} e^{i(m-m')\varphi} d\varphi = 2\pi \delta_{mm'} \quad (10) $$

$$ \int_0^{2\pi} \sin \varphi e^{i(m-m')\varphi} d\varphi = i\pi(\delta_{m,m'+1} - \delta_{m,m'-1}) \quad (11) $$

$$ \int_0^{2\pi} \cos \varphi e^{i(m-m')\varphi} d\varphi = \pi(\delta_{m,m'+1} + \delta_{m,m'-1}) \quad (12) $$

to obtain:

$$ \begin{aligned} & RQ \exp(-Qs^2 R^2 \sin^2 \theta - iR \cos \theta) \sin \theta (\delta_{m,1} + \delta_{m,-1}) \\ &= 2 \sum_{n=1}^{\infty} (-i)^{n+1} (2n+1) g_{n,TM}^m j_n(R) P_n^{|m|}(\cos \theta) \end{aligned} \quad (13) $$

$$ -iRQ \exp(-Qs^2 R^2 \sin^2 \theta - iR \cos \theta) \sin \theta (\delta_{m,1} - \delta_{m,-1}) \quad (14) $$

$$ = 2 \sum_{n=1}^{\infty} (-i)^{n+1} (2n+1) g_{n,TE}^m j_n(R) P_n^{|m|}(\cos \theta) $$

These Eqs. (13)–(14) imply that the non-zero values for the BSCs are for $m = \pm 1$, and that they satisfy a relation allowing one to define uni-index BSCs $g_n$ reading as:

$$ \frac{g_n}{2} = g_{n,TM}^1 = g_{n,TM}^{-1} = ig_{n,TE}^1 = -ig_{n,TE}^{-1} \quad (15) $$

This relation is typical of on-axis Gaussian beams, e.g. Eq. (6.3) in [1], more generally of on-axis axisymmetric beams of the first kind. Such beams are defined by the fact that the z-component $S_z$ of the time-averaged Poynting vector does not depend on the azimuthal coordinate $\varphi$, in a coordinate system in which the axis z is the symmetry axis of the beam. Such beams may be viewed as z-axisymmetric, e.g. Eq. (1) in [61], and Eqs. (6) and (7) in [62]. It has recently been demonstrated that such beams are actually circularly symmetric, i.e. all components of the Poynting vector does not depend on $\varphi$ [63].

As a result, Eqs. (13)–(14) may be expressed using a single equation reading as:

$$ RQ \exp(-Qs^2 R^2 \sin^2 \theta - iR \cos \theta) \sin \theta = \sum_{n=1}^{\infty} (-i)^{n+1} (2n+1) g_n j_n(R) P_n^1(\cos \theta) \quad (16) $$

In the next step, which is typical of the finite series approach, we specify a xy-plane location, i.e. $\theta = \pi/2$, leading to:

$$ R \exp(-s^2 R^2) = \sum_{n=1}^{\infty} (-i)^{n+1} (2n+1) g_n j_n(R) P_n^1(0) \quad (17) $$

When $\theta = \pi/2$, $P_n^m(\cos \theta)$ is nonzero only for $(n-m)$ even, which again is a typical result of the finite series approach in which the case $(n-m)$ odd requires another treatment (see Appendix A). For $(n-m)$ even, we then set $n = 2p+1$, leading to:

$$ x \exp(-s^2 x^2) = \sum_{p=0}^{\infty} a_p j_{2p+1}(x) \quad (18) $$
in which we conveniently change R (related to the physical problem in hand) to x (to insist on the mathematical problem in hand), and to:

$$a_p = (-1)^{p+1} (4p + 3) g_{2p+1} P_{2p+1}^1(0) \quad (19)$$

We then use [60]:

$$P_n^m(0) = \frac{(-1)^{(n+m)/2} (n+m-1)!!}{2^{(n-m)/2} (\frac{n-m}{2})!}, \quad (n-m) \text{ even} \quad (20)$$

and a classical expression for $n!! = 1.3.5...n$ reading as:

$$\Gamma(\frac{1}{2} + n) = \frac{(2n-1)!!}{2^n} \sqrt{\pi} \quad (21)$$

to obtain:

$$P_{2p+1}^{1}(0) = (-1)^{p+1} \frac{2\Gamma(p + 3/2)}{\sqrt{\pi} p!} \quad (22)$$

so that Eq. (19) becomes:

$$a_p = (4p + 3) \frac{2\Gamma(p + 3/2)}{\sqrt{\pi} p!} g_{2p+1} \quad (23)$$

The problem is then to evaluate the expansion coefficients $a_p$ of Eq. (18) from which we can obtain the uni-index BSCs $g_n$, $n$ odd, using Eq. (23).

## 3. First method

We provide two approaches to illustrate the first method used to evaluate the expansion coefficients.

### 3.1. First approach

Let us introduce the function $h(x)$ reading as:

$$h(x) = x^{3/2} \exp(-s^2 x^2) \quad (24)$$

and introduce its Maclaurin expansion reading as:

$$x^{3/2} \exp(-s^2 x^2) = \sum_{q=0}^{\infty} \frac{(-s^2)^q}{q!} x^{2q+3/2} \quad (25)$$

We next use a Bessel expansion reading as [64]:

$$x^q = \frac{2^{q+1}}{\sqrt{\pi}} \sum_{p=0}^{\infty} \frac{(q+2p+1/2)\Gamma(q+p+1/2)}{p!} j_{q+2p}(x) \quad (26)$$

leading to:

$$(x/2)^\mu = \frac{2}{\sqrt{\pi}} \sum_{k=0}^{\infty} \frac{(\mu + 2k + 1/2) \Gamma(\mu + k + 1/2)}{k!} j_{\mu + 2k}(x) \quad (27)$$

We may afterward use:

$$j_n(x) = \sqrt{\frac{\pi}{2x}} J_{n+1/2}(x) \quad (28)$$

in which $J_{n+1/2}(\cdot)$ denotes Bessel functions [65], allowing one to rewrite Eq. (27) as:

$$(x/2)^\mu = \sum_{k=0}^{\infty} \frac{(\mu + 2k)\Gamma(\mu + k)}{k!} J_{\mu + 2k}(x) \quad (29)$$

We may then return to Eq. (25) and express $x^{2q+3/2}$ using the Bessel series of Eq. (29) to obtain:

$$x^{3/2} \exp(-s^2 x^2) = 2^{3/2} \sum_{q=0}^{\infty} \sum_{k=0}^{\infty} \frac{(2q + 2k + 3/2)\Gamma(2q + k + 3/2)}{k!} \frac{(-4s^2)^q}{q!} \times J_{2q + 2k + 3/2}(x) \quad (30)$$

which, with $p = q + k$, may be rewritten as:

$$\begin{aligned} & x^{3/2} \exp(-s^2 x^2) \\ &= 2^{3/2} \sum_{p=0}^{\infty} \sum_{q=0}^{p} \frac{(2p + 3/2) \Gamma(q + p + 3/2)}{(p - q)!} \frac{(-4s^2)^q}{q!} J_{2p + 3/2}(x) \end{aligned} \quad (31)$$

$$= 2^{3/2} \sum_{p=0}^{\infty} (2p + 3/2) \sum_{q=0}^{p} \frac{\Gamma(p + q + 3/2)}{q!(p - q)!} (-4s^2)^q J_{2p+3/2}(x)$$

Using Eq. (28), Eq. (31) may be given the form of Eq. (18) with:

$$a_p = \frac{4}{\sqrt{\pi}} (2p + 3/2) \sum_{q=0}^{p} \frac{\Gamma(p + q + 3/2)}{q!(p - q)!} (-4s^2)^q \quad (32)$$

which is the first solution to the evaluation of the expansion coefficients $a_p$. The correctness of this result has been numerically checked with high precision calculations as follows. Relying on Eq. (18), and defining:

$$f(x) = \left| x \exp(-s^2 x^2) - \sum_{p=0}^{p_{\max}} a_p J_{2p+1}(x) \right| \quad (33)$$

it is indeed found that $f(x) \to 0$ when $p_{\max} \to \infty$ for various values of $s$ and $x$.

### 3.2. Second approach

Let us again consider Eq. (18) and implement Eq. (28) to obtain:

$$x^{3/2} \exp(-s^2 x^2) = \sqrt{\frac{\pi}{2}} \sum_{p=0}^{\infty} a_p J_{2p+3/2}(x) \quad (34)$$

We are now going to apply the Neumann Expansion Theorem (NET) which is the basic theorem used in the finite series technique to the evaluation of BSCs, see Chapter 5 in [1] and [64]. For this, Eq. (34) must be given the form of Eq. (5.2) in [1] reading as:

$$x^{1/2} g(x) = \sum_{n=0}^{\infty} c_n J_{n+1/2}(x) \quad (35)$$

This is done by introducing:

$$g(x) = x \exp(-s^2 x^2) \quad (36)$$

so that Eq. (34) becomes:

$$x^{1/2} g(x) = \sqrt{\frac{\pi}{2}} \sum_{p=0}^{\infty} a_p J_{2p+3/2}(x) \quad (37)$$

which is to be compared with Eq. (35). Setting $n = 2p + 1$ in the r.h.s. of Eq. (37) and equating the r.h.s. of Eq. (35), we obtain:

$$\sqrt{\frac{\pi}{2}} \sum_{n=1, \text{odd}}^{\infty} a_{(n-1)/2} J_{n+1/2}(x) = \sum_{n=0}^{\infty} c_n J_{n+1/2}(x) \quad (38)$$

so that:

$$\left\{ \begin{array}{l} c_n = 0, \text{ n even} \\ c_n = \sqrt{\frac{\pi}{2}} a_{(n-1)/2}, \text{ n odd} \end{array} \right. \quad (39)$$

that is to say:

$$\left\{ \begin{array}{l} c_{2p} = 0 \\ c_{2p+1} = \sqrt{\frac{\pi}{2}} a_p \end{array} \right. \quad (40)$$

The next step requires us to obtain the expansion coefficients $b_n$ defined by the expansion, see Eq. (5.3) in [1]:

$$g(x) = \sum_{n=0}^{\infty} b_n x^n \quad (41)$$

From Eqs. (36) and (41), we then have:

$$g(x) = \sum_{n=0}^{\infty} \frac{(-s^2)^n x^{2n+1}}{n!} = \sum_{n=0}^{\infty} b_n x^n \quad (42)$$

leading to:

$$\left\{ \begin{array}{l} b_n = 0, \text{ n even} \\ b_n = (-s^2)^{(n-1)/2} / [(n-1)/2]!, \text{ n odd} \end{array} \right. \quad (43)$$
The relationship between the coefficients $c_n$ of Eqs. (39)–(40) and the coefficients $b_n$ of Eq. (43) are given by the NET, see Eq. (5.4) of [1]:

$$c_n = (n + \frac{1}{2}) \sum_{m=0}^{\le n/2} 2^{n-2m+1/2} \frac{\Gamma(n-m+1/2)}{m!} b_{n-2m} \quad (44)$$

Since we only deal with non-zero values of $c_n$, i.e. those for $n$ odd as shown in Eq. (40), we have:

$$c_{2p+1} = \sqrt{\frac{\pi}{2}} a_p = (2p + \frac{3}{2}) \sum_{m=0}^{p} 2^{2p-2m+3/2} \frac{\Gamma(2p - m + 3/2)}{m!} b_{2p-2m+1} \quad (45)$$

which, setting $q = (p - m)$ becomes:

$$c_{2p+1} = \sqrt{\frac{\pi}{2}} a_p = (2p + \frac{3}{2}) \sum_{q=0}^{p} 2^{2q+3/2} \frac{\Gamma(p+q+3/2)}{(p-q)!} b_{2q+1} \quad (46)$$

in which we only need the coefficients $b_n$ with $n$ odd, rewritten from Eq. (43) with $n = 2q + 1$, reading as:

$$b_{2q+1} = (-s^2)^q / q! \quad (47)$$

Inserting Eq. (47) into Eq. (46), we obtain:

$$\sqrt{\frac{\pi}{2}} a_p = 2^{3/2} (2p + \frac{3}{2}) \sum_{q=0}^{p} \frac{\Gamma(p+q+3/2)(-4s^2)^q}{q!(p-q)!} \quad (48)$$

which identifies with Eq. (32) as it should.

## 4. Second method

We multiply Eq. (18) by $j_{2n+1}(x)$ and integrate from 0 to $\infty$, leading to:

$$\int_{0}^{\infty} x \exp(-s^2 x^2) j_{2n+1}(x) dx = \int_{0}^{\infty} \sum_{p=0}^{\infty} a_p j_{2p+1}(x) j_{2n+1}(x) dx \quad (49)$$

Interchanging the integral and the summation, we obtain:

$$\int_{0}^{\infty} x \exp(-s^2 x^2) j_{2n+1}(x) dx = \sum_{p=0}^{\infty} a_p \int_{0}^{\infty} j_{2p+1}(x) j_{2n+1}(x) dx \quad (50)$$

so that, from Eqs. (49) and (50), we obtain:

$$\int_{0}^{\infty} \sum_{p=0}^{\infty} a_p j_{2p+1}(x) j_{2n+1}(x) dx = \sum_{p=0}^{\infty} a_p \int_{0}^{\infty} j_{2p+1}(x) j_{2n+1}(x) dx \quad (51)$$

But we have [65]:

$$\int_{0}^{\infty} j_{n}(x) j_{m}(x) dx = \frac{\sin[(n - m)(\pi/2)]}{n(n + 1) - m(n + 1)} \quad \text{for } n \neq m \quad (52)$$

$$\int_{0}^{\infty} j_{n}(x) j_{m}(x) dx = \frac{\pi}{2(2n + 1)} \quad \text{for } n = m \quad (53)$$

In Eq. (51), we have to deal with $n = 2p + 1$ and $m = 2n + 1$, so that Eqs. (52) and (53) reduce to:

$$\int_{0}^{\infty} j_{2p+1}(x) j_{2n+1}(x) dx = \frac{\pi}{2(4n + 3)} \delta_{np} \quad (54)$$

Inserting Eq. (54) into Eq. (50), we readily obtain:

$$a_p = \frac{2(4p + 3)}{\pi} \int_{0}^{\infty} x \exp(-s^2 x^2) j_{2p+1}(x) dx \quad (55)$$

which is the second solution to the evaluation of the expansion coefficients $a_p$.

## 5. Discussion

The two solution of Eqs. (32) and (55) are compared in the figure, using logarithmic scales. Furthermore, absolute values are used because the series may become alternating, e.g. the Table in Appendix C. The results for four values of the beam confinement factor $s$ are displayed, namely 1/10, 1/13, 1/17 and 1/20. We observe that the first solution

of Eq. (32) displayed with a blue line and blue symbols blows up in contrast with the second solution of Eq. (55). These behaviors may be checked analytically, see Appendix B.

The ranges chosen for $p$ depend on $s$ and have been chosen to visually exhibit the blowing-ups, from 0 to 40 when $s = 1/10$, up to from 0 to more than 80 for $s = 1/20$. From [54], we extract two successful criteria to the evaluation of the critical values $n_c = 2p_c + 1$ for the onset of the blowing-ups (see Fig. 1).

One criterion is theoretical and reads as $p_c = \text{int}[e/(16s^2)]$ leading to $p_c = 17, 29, 49$ and $68$ for $s = 1/10, 1/13, 1/17$ and $1/20$ respectively. The second criterion is empirical and can read $p_c = \text{int}[3/(16s^2) - 1/2]$ leading to $p_c = 18, 31, 54$ and $74$ respectively. Both criteria agree reasonably well with the critical values available from the figure. It should be noted that, for $p < p_c$, the two solutions seem to provide the same results. This is actually not true. The results are different although the differences are so small that they cannot be seen on the figure, see however the Table in Appendix C. These differences are 0 only for $s = 0$.

From a physical point of view the blowing-ups correspond to the presence of evanescent waves, see [54,55], [56], while the other results correspond to the elimination of the evanescent waves (at least most of them).

The reason for the disagreement between the two solutions is the fact that the interchange between the integral and the summation carried out in the step from Eq. (49) is not allowed (see Appendix C). By nevertheless carrying out the interchange, we therefore force the mathematics to produce a solution in which the interchange is allowed, and “indirectly” eliminate the evanescent waves (or at least most of them). We may demonstrate that this last interchange is indeed allowed using the Dominated Convergence Theorem [66], but such a demonstration may be omitted since it is a direct consequence of the procedure used to obtain the second solution of Eq. (55), in which Eq. (51) is indeed used to force it.

As a last word, let us mention that, although the obtained second solution is physically of interest, we have to insist on the fact that it is mathematically incorrect. Indeed, we may numerically check that Eq. (18) is not satisfied when we use Eq. (55) to specify the values of $a_p$'s. This is maybe the most important message of the present paper, namely that a mathematically incorrect procedure may be used to produce a physically correct solution of interest. Actually, any set of BSCs generates a genuine electromagnetic field, even if these coefficients were tossed, and therefore generates a physically correct solution. But, in the present case, the solution obtained is not only physically correct but it is also of interest since it eliminates the evanescent waves (or at least most of them).

## 6. Conclusion

Usually, physicists do not care too much with interchanges between integrals and summations. This is usually justified by the fact that equations, in physics, most often do not exhibit “pathological” features which would prevent such interchanges to be carried out. The present work provides a counter-example.

But there is an interesting consequence to be noted. Namely, by using an “incorrect” mathematical step, we force the mathematics to provide a solution which is physically a “correct” solution, by eliminating evanescent waves (or at least most of them).

It must be furthermore noted that the blowing-ups observed with the correct solution do not imply a blowing-up of the associated physical quantities, in particular do not imply a blowing-up of the associated fields. This is because, actually, when dealing with physical quantities, BSCs which blow-up are multiplied by functions decreasing to zero fast enough to compensate for the blowing-ups. We then have obtained two solutions which, both of them, are physically correct. Comparisons between physical quantities, using either of the two solutions would then allow one to characterize the influence of the evanescent waves [67].
<figure><img src="image_4.png" /><figcaption>Fig. 1. Expansion coefficients $a_p$ for various values of the beam confinement parameter. Blue: blowing-up solution. Black: blowing-up eliminated.</figcaption></figure>

## CRediT authorship contribution statement

**Gérard Gouesbet:** Writing – original draft, Formal analysis, Data curation, Conceptualization. **Jianqi Shen:** Writing – review & editing, Formal analysis, Data curation, Conceptualization. **Leonardo André Ambrosio:** Writing – review & editing, Funding acquisition, Formal analysis, Data curation, Conceptualization.

## Funding information

L. Ambrosio thanks the Council for Scientific and Technological Development (CNPq), Brazil (309201/2021-7) and the São Paulo Research Foundation (FAPESP), Brazil (2020/05280-5) for partially supporting this work.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

We are pleased to thank Bernard Julia, from the “École Normale Supérieure in Paris, France”, Richard Panetta, from “Texas A&M University, USA”, and André Draux from “Rouen INSA, France” for advices which helped us to achieve the present work.

## Appendix A

The paper deals with the case (n-m) even. The case (n-m) odd is not considered in this paper because it would behave in a completely similar way, and then would unnecessarily increase the size of the paper. The similarity between the two cases is obvious from two equations below, providing the expressions of BSCs using finite series

for on-axis Gaussian beams at an order of approximation called L⁻, see [1], page 167:

$$g_{2p+1} = \sum_{j=0}^{p} \frac{p!}{j!(p-j)!} \frac{\Gamma(p + j + 3/2)}{\Gamma(p + 3/2)} (-4s^2)^j \quad (56)$$

$$g_{2p+2} = \sum_{j=0}^{p} \frac{p!}{j!(p-j)!} \frac{\Gamma(p + j + 5/2)}{\Gamma(p + 5/2)} (-4s^2)^j \quad (57)$$

## Appendix B

In a first step, we analytically demonstrate that the second solution of Eq. (55) does not blow up. Using Eq. (28), Eq. (55) may be rewritten as:

$$a_p = (4p + 3) \frac{\sqrt{2}}{\pi} \int_0^\infty \exp(-s^2 x^2) J_{p+3/2}(x) x^{1/2} dx \quad (58)$$

From Eq. (6.631.1) of [68], with Re(a) > 0 and Re(n + t) > -1, we have:

$$\begin{aligned} & \int_{0}^{\infty} \exp(-ax^2) J_n(bx) x' dx \\ &= \frac{b^n \Gamma(\frac{n+t+1}{2})}{2^{n+1} a^{(n+t+1)/2} \Gamma(n+1)} [I_1(\frac{n+t+1}{2}; n+1; \frac{-b^2}{4a})] \\ &= \frac{b^n \Gamma(\frac{n+t+1}{2})}{2^{n+1} a^{(n+t+1)/2} \Gamma(n+1)} \exp(\frac{-b^2}{4a}) I_1(\frac{n-t+1}{2}; n+1; \frac{b^2}{4a}) \end{aligned} \quad (59)$$

in which $I_1$ is a hypergeometric function, also called the Kummer's function of the first kind, and in which we have used the Kummer's transformation [69]:

$$I_1 F_1(a; b; x) = \exp(x) I_1 F_1(b - a; b; -x) \quad (60)$$

Letting $a = s^2$, $b = 1$, $n = 2p + 3/2$ and $t = 1/2$, we then obtain:

$$a_p = \frac{2^{-2p-1} s^{-2p-3}}{\sqrt{\pi}} \frac{\Gamma(p + 3/2)}{\Gamma(2p + 3/2)} \exp\left(\frac{-1}{4s^2}\right) I_1\left(p + 1; 2p + 5/2; \frac{1}{4s^2}\right) \quad (61)$$
We then use a classical expression reading as:

$$ {}_1F_1(a; b; z) = \sum_{q=0}^{\infty} \frac{(a)_q z^q}{(b)_q q!} \qquad (62) $$

in which $(c)_q$ is the Pochhammer symbol defined as:

$$ (c)_0 = 1 \qquad (63) $$

$$ (c)_q = c(c+1)(c+2)...(c+q-1), \quad q \neq 0 \qquad (64) $$

so that the hypergeometric function of Eq. (61) may be rewritten as:

$$ {}_1F_1(p+1; 2p+5/2; \frac{1}{4s^2}) = \sum_{q=0}^{\infty} \frac{(p+1)_q}{(2p+5/2)_q} \left(\frac{1}{4s^2}\right)^q \frac{1}{q!} \qquad (65) $$

Similarly:

$$ {}_1F_1(p+1; 2p+2; \frac{1}{4s^2}) = \sum_{q=0}^{\infty} \frac{(p+1)_q}{(2p+2)_q} \left(\frac{1}{4s^2}\right)^q \frac{1}{q!} \qquad (66) $$

from which we readily deduce:

$$ {}_1F_1(p+1; 2p+5/2; \frac{1}{4s^2}) < {}_1F_1(p+1; 2p+2; \frac{1}{4s^2}) \qquad (67) $$

But, from [69], we have:

$$ {}_1F_1(v, 2v, z) = \Gamma(v+1/2) \exp(z/2)(z/4)^{-v+1/2} I_{v-1/2}(z/2) \qquad (68) $$

in which $I_\alpha(x)$ is the modified Bessel function satisfying the relation:

$$ I_{\alpha}(x) = \sum_{s=0}^{\infty} \frac{1}{s!(s+\alpha)!} \left(\frac{x}{\alpha}\right)^{2s+\alpha} \qquad (69) $$

Using Eqs. (61) and (68), and the inequality of (67), we may then establish another inequality reading as:

$$ a_p < b_p = \frac{2^{2p+1} s^{-2} \left[ \frac{\Gamma(p+3/2)}{\Gamma(2p+3/2)} \right]^2}{\sqrt{\pi} \Gamma(2p+3/2)} \exp\left(\frac{-1}{8s^2}\right) I_{p+1/2}\left(\frac{1}{8s^2}\right) \qquad (70) $$

Using $\Gamma(x+1) = x\Gamma(x)$ and Eq. (69), we finally find that, for $p$ sufficiently large:

$$ \frac{b_{p+1}}{b_p} \approx \frac{I_{p+3/2}(\frac{1}{8s^2})}{I_{p+1/2}(\frac{1}{8s^2})} < 1 \qquad (71) $$

meaning that the $b_p$s, and therefore the $a_p$s according to Eq. (70), do not blow up.

In a second step, we analytically demonstrate that the first solution of Eq. (32) indeed blows up. To begin with, this equation is rewritten as:

$$ \begin{aligned} a_p &= \frac{2(4p+3)}{\sqrt{\pi}} \frac{\Gamma(p+3/2)}{p!} \sum_{q=0}^{p} \binom{p}{q} \frac{\Gamma(p+q+3/2)}{\Gamma(p+3/2)} (-4s^2)^q \\ &= \frac{2(4p+3)}{\sqrt{\pi}} \frac{\Gamma(p+3/2)}{p!} \sum_{q=0}^{p} \binom{p}{q} (p+3/2)_q (-4s^2)^q \\ &= \frac{4}{p! \sqrt{\pi}} \frac{4p+3}{2p+3} \Gamma(p+5/2) \sum_{q=0}^{p} \binom{p}{q} (p+3/2)_q (-4s^2)^q \end{aligned} \qquad (72) $$

in which $\binom{p}{q}$ is the binomial coefficient and $(a)_q$ is again the Pochhammer coefficient.

We then introduce the generalized Bessel polynomial $y_p(x; a)$ which, according to Grosswald [70], see as well Eq. (29) in [55], reads as:

$$ y_p(x; a) = \sum_{q=0}^{p} \binom{p}{q} (p+a-1)_q \left(\frac{x}{2}\right)^q \qquad (73) $$

allowing one to deal again with Eq. (72) to rewrite it as:

$$ \begin{aligned} a_p &= \frac{2(4p+3)}{\sqrt{\pi}} \frac{\Gamma(p+3/2)}{p!} \sum_{q=0}^{p} \binom{p}{q} (p+3/2)_q (-4s^2)^q \\ &= \frac{2(4p+3)}{\sqrt{\pi}} \frac{\Gamma(p+3/2)}{p!} y_p(-8s^2; 5/2) \end{aligned} \qquad (74) $$

But, when $p$ is large, $y_p(x; a)$ tends to the expression of Eq. (75) given below, see again Grosswald [70] and as well Eq. (31) in [55]:

$$ y_p(x; a) = \left(\frac{2x}{e} p\right)^p 2^{a-3/2} e^{1/x} [1 + O(1/p)] \qquad (75) $$

so that, for $p$ large enough, we have:

$$ a_p \approx \frac{4(4p+3)}{\sqrt{\pi}} \frac{\Gamma(p+3/2)}{p!} \left(\frac{-16s^2}{e} p\right)^p \exp\left(\frac{-1}{8s^2}\right) \qquad (76) $$

leading to:

$$ \frac{a_{p+1}}{a_p} \approx \frac{-16s^2}{e} p \qquad (77) $$

implying the blowing-up of the coefficients $a_p$.

## Appendix C

Let $IS$ be:

$$ IS = \int_{0}^{\infty} \sum_{p=0}^{\infty} a_p j_{2p+1}(R) j_{2n+1}(R) dR \qquad (78) $$

and $SI$ be:

$$ SI = \sum_{p=0}^{\infty} a_p \int_{0}^{\infty} j_{2p+1}(R) j_{2n+1}(R) dR \qquad (79) $$

both quantities being taken from Eq. (51). In $IS$, the summation $S$ is carried out first and the integral $I$ afterward. In $SI$, the integral is carried out first and the summation afterward.

We also define the difference:

$$ D = IS - SI \qquad (80) $$

All these quantities are displayed in the Table below, using for $a_p$ the first solution of Eq. (32), setting the beam confinement parameter $s$ to $1/10$, and running the computations with 3000 digits.

<table><thead><tr><th>p</th><th>IS</th><th>SI</th><th>D</th></tr></thead><tbody><tr><td>1</td><td>2.120575</td><td>2.120575</td><td>1.38E-10</td></tr><tr><td>2</td><td>2.194795</td><td>2.194795</td><td>-1.65E-10</td></tr><tr><td>3</td><td>1.953446</td><td>1.953446</td><td>2.11E-10</td></tr><tr><td>4</td><td>1.547777</td><td>1.547777</td><td>-2.93E-10</td></tr><tr><td>5</td><td>1.109633</td><td>1.109633</td><td>4.36E-10</td></tr><tr><td>6</td><td>0.72692</td><td>0.72692</td><td>-6.97E-10</td></tr><tr><td>7</td><td>0.438198</td><td>0.438198</td><td>1.20E-09</td></tr><tr><td>8</td><td>0.244416</td><td>0.244416</td><td>-2.19E-09</td></tr><tr><td>9</td><td>0.126731</td><td>0.126731</td><td>4.29E-09</td></tr><tr><td>10</td><td>0.061333</td><td>0.061333</td><td>-8.94E-09</td></tr><tr><td>11</td><td>0.027806</td><td>0.027806</td><td>1.98E-08</td></tr><tr><td>12</td><td>0.011848</td><td>0.011848</td><td>-4.65E-08</td></tr><tr><td>13</td><td>0.004759</td><td>0.004759</td><td>1.15E-07</td></tr><tr><td>14</td><td>0.001807</td><td>0.001807</td><td>-3.02E-07</td></tr><tr><td>15</td><td>0.00065</td><td>0.000649</td><td>-8.34E-07</td></tr><tr><td>16</td><td>0.000222</td><td>0.000225</td><td>-2.42E-06</td></tr><tr><td>17</td><td>7.23E-05</td><td>6.49E-05</td><td>7.37E-06</td></tr><tr><td>18</td><td>2.24E-05</td><td>4.59E-05</td><td>-2.349E-05</td></tr><tr><td>19</td><td>6.66E-06</td><td>-7.2E-05</td><td>7.8332E-05</td></tr><tr><td>20</td><td>1.89E-06</td><td>0.000275</td><td>-0.0002727</td></tr><tr><td>21</td><td>5.15E-07</td><td>-0.00099</td><td>0.00098976</td></tr><tr><td>22</td><td>1.35E-07</td><td>0.00374</td><td>-0.0037395</td></tr><tr><td>23</td><td>3.39E-08</td><td>-0.01469</td><td>0.01468854</td></tr><tr><td>24</td><td>8.20E-09</td><td>0.059904</td><td>-0.0599044</td></tr><tr><td>25</td><td>1.29E-09</td><td>-0.25336</td><td>0.25335919</td></tr><tr><td>26</td><td>4.32E-10</td><td>1.109987</td><td>-1.1099872</td></tr><tr><td>27</td><td>9.41E-11</td><td>-5.03193</td><td>5.031922733</td></tr><tr><td>28</td><td>1.98E-11</td><td>23.5799</td><td>-23.579905</td></tr><tr><td>29</td><td>4.05E-12</td><td>-114.109</td><td>114.109102</td></tr></tbody></table>

We observe that, for $p > 1$, $IS$ decreases steadily when $p$ increases while $SI$ decreases steadily up to $p = 18$. From this value $p = 18$, $SI$
becomes alternating with absolute values increasing steadily, exhibiting a blowing-up. For small values of $p$ up to 14, IS and SI perfectly agree up to 6 digits and afterward disagree more and more. We recall that we found that the critical value $p_c$ was found to be 17 using the theoretical criterion and 18 using the empirical criterion. These values are well confirmed by the table. Apart of these remarks, the most important message is that interchanging the integral and the summation as in Eq. (51) is indeed forbidden. Doing this, we then generate another solution, which eliminates the evanescent waves, and for which the interchange is allowed since it is forced to obtain the second solution, which does not blow-up.

## Data availability

No data was used for the research described in the article.

## References

[1] Gouesbet G, Gréhan G. Generalized Lorenz-Mie theories. 3rd ed.. Springer; 2023.

[2] Gouesbet G, Maheu B, Gréhan G. Light scattering from a sphere arbitrarily located in a Gaussian beam, using a Bromwich formulation. J Opt Soc Amer A 1988;5:9-1427-43.

[3] Maheu B, Gouesbet G, Gréhan G. A concise presentation of the generalized Lorenz-Mie theory for arbitrary location of the scatterer in an arbitrary incident profile. J Optics (Paris) 1988;19:2-59-67.

[4] Waterman PC. Symmetry, unitarity, and geometry in electromagnetic scattering. Phys Rev D 1971;3: 4:825-39.

[5] Mishchenko MI, Travis LD, Lacis AA. Scattering, absorption, and emission of light by small particles. Cambridge: Cambridge University Press; 2002.

[6] Mishchenko MI. Electromagnetic scattering by particles and particle groups, an introduction. Cambridge, UK: Cambridge University Press; 2014.

[7] Mackowski DW, Mishchenko MI. Direct simulation of multiple scattering by discrete random media illuminated by Gaussian beams. Phys Rev A 2011;83:013804.

[8] Wang J, Chen A, Han Y, Briard P. Light scattering from an optically anisotropic particle illuminated by an arbitrary shaped beam. J Quant Spectrosc Radiat Transfer 2015;167:135-44.

[9] Gouesbet G. Latest achievements in generalized Lorenz-Mie theories: A commented reference database. Ann Phys, Lpz 2014;526, 11-12:461-89.

[10] Gouesbet G. T-matrix methods for electromagnetic structured beams: A commented reference database for the period 2014-2018. J Quant Spectrosc Radiat Transfer 2019;230:247-81.

[11] Gouesbet G. T-matrix methods for electromagnetic structured beams: A commented reference database for the period 2019-2023. J Quant Spectrosc Radiat Transfer 2024;322:109015.

[12] Gouesbet G, Letellier C, Ren KF, Gréhan G. Discussion of two quadrature methods of evaluating beam shape coefficients in generalized Lorenz-Mie theory. Appl Opt 1996;35:9-1537-42.

[13] Gouesbet G, Ambrosio LA, Lock JA. On an infinite number of quadratures to evaluate beam shape coefficients in generalized Lorenz-Mie theory and extended boundary condition method for structured EM fields. J Quant Spectrosc Radiat Transfer 2020;242:196779.

[14] Gouesbet G, Lock JA, Gréhan G. Generalized Lorenz-Mie theories and description of electromagnetic arbitrary shaped beams: localized approximations and localized beam models, a review. J Quant Spectrosc Radiat Transfer 2011;112:1-27.

[15] Wang J, Gouesbet G. Note on the use of localized beam models for light scattering theories in spherical coordinates. Appl Opt 2012;51, 17:3832-6.

[16] Gouesbet G. Second modified localized approximation for use in generalized Lorenz-Mie theories and other theories revisited. J Opt Soc Amer A 2013;30: 4:560-4.

[17] Gouesbet G, Lock JA. Comments on localized and integral localized approximations in spherical coordinates. J Quant Spectrosc Radiat Transfer 2016;179:132-6.

[18] Gouesbet G. On the validity of localized approximations for Bessel beams: All N-Bessel beams are identically equal to zero. J Quant Spectrosc Radiat Transfer 2016;176:82-6.

[19] Gouesbet G, Lock JA, Ambrosio LA, Wang J. On the validity of localized approximation for an on-axis zeroth-order Bessel beam. J Quant Spectrosc Radiat Transfer 2017;195:18-25.

[20] Ambrosio LA, Wang J, Gouesbet G. On the validity of the integral localized approximation for Bessel beams and associated radiation pressure forces. Appl Opt 2017;56, 19:5377-87.

[21] Chafiq A, Ambrosio LA, Gouesbet G, Belafhal A. On the validity of the integral localized approximation for on-axis zeroth-order Mathieu beams. J Quant Spectrosc Radiat Transfer 2018;204:27-34.

[22] Gouesbet G, Ambrosio LA. On the validity of the use of a localized approximation for helical beams. I. Formal aspects. J Quant Spectrosc Radiat Transfer 2018;208:12-8.

[23] Ambrosio LA, Gouesbet G. On the validity of the use of a localized approximation for helical beams. II. Numerical aspects. J Quant Spectrosc Radiat Transfer 2018;215:41-50.

[24] Ambrosio LA, Votto LFM, Gouesbet G, Wang J. Assessing the validity of the localized approximation for discrete superposition of Bessel beams. J Opt Soc Am B 2018;35, 11:2690-8.

[25] Ambrosio LA, Gouesbet G. On localized approximations for Laguerre-Gauss beams focused by a lens. J Quant Spectrosc Radiat Transfer 2018;218:100-14.

[26] Gouesbet G, Lock JA, Han Y, Wang J. Efficient computation of arbitrary beam scattering on a sphere: Comments and rebuttal, with a review on the angular spectrum decomposition. J Quant Spectrosc Radiat Transfer 2021;276:107913.

[27] Qiu J, Shen J. Beam shape coefficient calculation for a Gaussian beam: localized approximation, quadrature and angular spectrum decomposition methods. Appl Opt 2018;57, 2:302-13.

[28] Shen J, Liu X, Wang W, Yu H. Calculation of light scattering of an elliptical Gaussian beam by a spherical particle. J Opt Soc Amer A 2018;35, 8:1288-98.

[29] Wang C, Shen J, Ren K-F. Spherical harmonics expansion of the evanescent waves in angular spectrum decomposition of shaped beams. J Quant Spectrosc Radiat Transfer 2020;251:107012.

[30] Cheng Y, Shen J. Curved photonic jet produced from a spherical particle illuminated by a Bessel-Gaussian beam. J Quant Spectrosc Radiat Transfer 2021;272:107765.

[31] Zhang Y, Shen J. Comparison of a standard elliptical Bessel beam and a refracted circular Bessel beam at oblique incidence. J Quant Spectrosc Radiat Transfer 2021;272:107773.

[32] Ren KF, Shen J. In: Kokhanovsky A, editor. Volume 7: Light absorption and scattering in turbid media. Springer; 2021, p. 125-58.

[33] Shen J, Liu J, Wang Y, Liu Z, Yu H. Cylindrical wave spectrum decomposition method for evaluating the expansion coefficients of the shaped beam in spherical coordinates. J Quant Spectrosc Radiat Transfer 2022;283:108138.

[34] Shen J, Wang Y, Yu H, Ambrosio LA, Gouesbet G. Angular spectrum representation of the Bessel-Gauss beam and its approximation: A comparison with the localized approximation. J Quant Spectrosc Radiat Transfer 2022;284:108167.

[35] Shen J, Liu J, Liu Z, Yu H. Angular spectrum decomposition method and the quadrature method in the generalized Lorenz-Mie theory for evaluating the beam shape coefficients of TEM01* doughnut beam. Opt Commun 2022;515:128224.

[36] Liu Z, Shen J, Yu H. Scattering of Laguerre-Gauss light beam by a sphere: The angular spectrum decomposition method and a comparison with the localized approximation method. J Quant Spectrosc Radiat Transfer 2022;287:108214.

[37] Shen J, Liu Z, Liu J. Expanding non-axisymmetric beams in spherical coordinates with cylindrical wave spectrum decomposition. Results Phys 2022. http://dx.doi.org/10.1016/j.rinp.2022.106068.

[38] Zhong S, Lin J, Shen J. Beam shape coefficients of hollow vortex Gaussian beam and near field scattering. J Opt Soc Amer A 2024;41(7):1403-12.

[39] Shen J, Yu H. Radial quadrature method for evaluating the beam shape coefficients in spherical coordinates. J Quant Spectrosc Radiat Transfer 2023;305:108627.

[40] Lin J, Zhong S, Shen J. Equivalence between radial quadrature and finite series for spherical wave expansion of Bessel beams. J Opt Soc Amer A 2023;40(6):1201-7.

[41] Tang S, Wang M, Shen J. Evaluating the beam shape coefficients of Bessel-Gauss beams with radial quadrature: a comparison with angular spectrum decomposition and finite series methods. J Opt Soc Amer A 2024;41(6):1152-65.

[42] Wang M, Tang S, Shen J. Radial quadrature method for evaluating the beam shape coefficients of Laguerre-Gauss beam. J Opt Soc Am B 2024;41(7):1587-98.

[43] Gouesbet G, Gréhan G, Maheu B. Expressions to compute the coefficients $g_n^m$ in the generalized Lorenz-Mie theory, using finite series. J Opt (Paris) 1988;19:1-35-48.

[44] Gouesbet G, Gréhan G, Maheu B. Computations of the $g_n^m$ coefficients in the generalized Lorenz-Mie theory using three different methods. Appl Opt 1988;27,23:4874-83.

[45] Gouesbet G, Votto LFM, Ambrosio LA. Finite series expressions to evaluate the beam shape coefficients of a Laguerre-Gauss beam freely propagating. J Quant Spectrosc Radiat Transfer 2019;227:12-9.

[46] Gouesbet G, Ambrosio LA, Votto LFM. Finite series expressions to evaluate the beam shape coefficients of a Laguerre-Gauss beam focused by a lens in an on-axis configuration. J Quant Spectrosc Radiat Transfer 2019;242:106759.

[47] Votto LFM, Ambrosio LA, Gouesbet G. Evaluation of beam shape coefficients of paraxial Laguerre-Gauss beam freely propagating by using three remodeling methods. J Quant Spectrosc Radiat Transfer 2019;239:106618.

[48] Ambrosio LA, Gouesbet G. Modified finite series technique for the evaluation of beam shape coefficients in the T-matrix methods for structured beams with application to Bessel beams. J Quant Spectrosc Radiat Transfer 2020;248:107007.

[49] Valdivia NL, Votto LFM, Gouesbet G, Wang J, Ambrosio LA. Bessel-Gauss beams in the generalized Lorenz-Mie theory using three remodeling techniques. J Quant Spectrosc Radiat Transfer 2020;256:107292.
[50] Votto LFM, Ambrosio LA, Gouesbet G, Wang J. Finite series algorithm design for lens-focused Laguerre-Gauss beams in the generalized Lorenz-Mie theory. J Quantitat Spectroscopy Radiat Transfer Virtual special issue of LIP2020 2021;261. Article 107 488.

[51] Votto LFM, Chafiq A, Belafhal A, Gouesbet G, Ambrosio LA. Hermite-Gaussian beams in the generalized Lorenz-Mie theory through finite-series Laguerre-Gaussian beam shape coefficients. J Opt Soc Am B 2022;39, 4:1027-32 with a supplemental document.

[52] Votto LFM, Chafiq A, Gouesbet G, Ambrosio LA, Belafhal A. Ince-Gaussian beams in the generalized Lorenz-Mie theory through finite series Laguerre-Gaussian beam shape coefficients. J Quant Spectrosc Radiat Transfer 2023;302:108565.

[53] Votto LFM, Gouesbet G, Ambrosio LA. A framework for the finite series method of the generalized Lorenz-Mie theory and its application to freely propagating Laguerre-Gaussian beams. J Quant Spectrosc Radiat Transfer 2023;309:108706.

[54] Gouesbet G, Votto LFM, Ambrosio LA. Evaluation of beam shape coefficients in T-matrix methods using a finite series technique: On blow-ups using hypergeometric functions and generalized Bessel polynomials. J Opt Soc Am B 2023;40(12):3174-81.

[55] Votto LFM, Gouesbet G, Ambrosio LA. Blowing-ups of beam shape coefficients of Gaussian beams using finite series in generalized Lorenz-Mie theory. J Quant Spectrosc Radiat Transfer 2023;311:108787.

[56] Shen J, Tang S, Ambrosio LA, Gouesbet G. On evanescent waves and blowing-ups of the finite series technique in spherical wave expansion of shaped beams. J Quant Spectrosc Radiat Transfer 2024;313:108846.

[57] Davis LW. Theory of electromagnetic beams. Phys Rev 1979;19, 3:1177-9.

[58] Gouesbet G, Lock JA, Gréhan G. Partial wave representations of laser beams for use in light scattering calculations. Appl Opt 1995;34,12:2133-43.

[59] Gouesbet G, Shen J, Ambrosio LA. Diverging and converging schemes of approximations to describe fundamental EM Gaussian beams beyond the paraxial approximation. J Quant Spectrosc Radiat Transfer 2022;291:108344.

[60] Robin L. Fonctions sphériques de Legendre et fonctions sphéroidales. volumes 1, 2, 3. Gauthier-Villars, Paris; 1957.

[61] Gouesbet G, Ambrosio LA. Rayleigh limit of the generalized Lorenz-Mie theory for on-axis beams and its relationship with the dipole theory of forces. Part i. non dark axisymmetric beams of the first kind, with the example of Gaussian beams. J Quant Spectrosc Radiat Transfer 2021;266:107569.

[62] Gouesbet G, Votto LFM, Ambrosio LA, Wang J. Poynting vector and beam shape coefficients: On new families of symmetries (non-dark axisymmetric beams of the second kind and dark axisymmetric beams). J Quant Spectrosc Radiat Transfer 2021;271:107745.

[63] Gouesbet G. Symmetries of circularly axisymmetric beams expressed in terms of beam shape coefficients. J Quant Spectrosc Radiat Transfer 2024;316:108919.

[64] Watson GN. A treatise of the theory of Bessel functions. 2nd ed.. Cambridge: Cambridge at the University Press; 1962.

[65] Arfken G. Mathematical methods for physicists. London, New-York: Academic Press; 1966, 1968.

[66] NIST handbook of Mathematical Functions.

[67] Tang S, Shen J, Gouesbet G, Ambrosio LA. On radial quadrature method applied to spherical wave expansion of Gaussian beams. J Quant Spect Radiat Trans 2024. [submitted for publication].

[68] Gradshteyn IS, Ryzhik IM. Table of integrals, series, and products. Academic Press; 1980.

[69] Polyanin AD, Manzhirov AV. Handbook of mathematics for engineers and scientists. Chapman and Hall; 2007.

[70] Grosswald E. Bessel polynomials, Lecture notes in mathematics. Springer-Verlag; 1978.