Electrostatic limit of the $T$-matrix for electromagnetic scattering: Exact results for spheroidal particles

TL;DR

This paper derives exact electrostatic expressions for the T-matrix of spheroidal particles, clarifying the long-wavelength limit of electromagnetic scattering and providing tools for analyzing convergence and the Rayleigh hypothesis.

Contribution

It introduces a fully electrostatic version of the T-matrix method (ES-EBCM), deriving analytic expressions for spheroids and connecting them to the electromagnetic T-matrix in the long-wavelength limit.

Findings

Explicit electrostatic T-matrix expressions for spheroids.

Symmetry properties of the ES-EBCM matrices.

Insights into the convergence of the T-matrix formalism.

Abstract

The $T$-matrix, often obtained with Waterman's extended boundary condition method (EBCM), is a widely-used tool for fast calculations of electromagnetic scattering by particles. Here we investigate the quasistatic or long-wavelength limit of this approach, where it reduces to an electrostatics problem. We first present a fully electrostatic version of the EBCM/T-matrix method (dubbed ES-EBCM). Explicit expressions are then given to quantitatively express the long-wavelength limit of the EBCM matrix elements in terms of those of the ES-EBCM formalism. From this connection we deduce a number of symmetry properties of the ES-EBCM matrices. We then investigate the matrix elements of the ES-EBCM formalism in the special case of prolate spheroids. Using the general electrostatic solution in spheroidal coordinates, we derive fully analytic expressions (in the form of finite sums) for all matrix elements. Those can be used for example for studies of the convergence of the $T$-matrix formalism. We illustrate this by discussing the validity of the Rayleigh hypothesis, where analytical expressions highlight clearly the link with analytical continuation of series.

Equations

$$\mathbf{E}(\mathbf{r})=\sum _{n,m}{a}_n^m{\mathbf{M}}_{nm}(\mathbf{r})+b_n^m{\mathbf{N}}_{nm}(\mathbf{r})$$

$$G(\mathbf{r},{\mathbf{r}}^{\prime })=$$

$${\int }_S[{\Phi }_{\mathrm{in}}({\mathbf{r}}^{\prime })\frac{\partial G}{\partial n_{-}^{\prime }}$$

$$=\{\begin{array}{ll}-{\Phi }_{\mathrm{in}}(\mathbf{r})& \text{mbox}if\mathbf{r}\in V,\\ 0& \text{mbox}if\mathbf{r}\notin V.\end{array}$$

$${\int }_S[{\Phi }_{\mathrm{out}}({\mathbf{r}}^{\prime })\frac{\partial G}{\partial n_{+}^{\prime }}$$

$$=\{\begin{array}{ll}{\Phi }_{\mathrm{sca}}(\mathbf{r})& \text{mbox}if\mathbf{r}\notin V,\\ -{\Phi }_{\mathrm{ext}}(\mathbf{r})& \text{mbox}if\mathbf{r}\in V.\end{array}$$

$$\left(ϵ-1\right){\int }_S\frac{\partial {\Phi }_{\mathrm{in}}}{\partial n_{-}^{\prime }}({\mathbf{r}}^{\prime })$$

$$=\{\begin{array}{ll}{\Phi }_{\mathrm{ext}}-{\Phi }_{\mathrm{in}},& \text{mbox}if\mathbf{r}\in V\\ -{\Phi }_{\mathrm{sca}},& \text{mbox}if\mathbf{r}\notin V\end{array}$$

$${\Psi }_{nm}^{(1)}(\mathbf{r})=$$

$${\Psi }_{nm}^{(3)}(\mathbf{r})=$$

$$A_n^m=$$

$${\int }_{r=R}{\Psi }_{nm}^{(i)}(\mathbf{r}){\Psi }_{n^{\prime }{m}^{\prime }}^{(j)}(\mathbf{r})dS=\frac{{\delta }_{n,n^{\prime }}{\delta }_{m,m^{\prime }}}{2n+1}{R}^2,$$

$$G(\mathbf{r},{\mathbf{r}}^{\prime })=\{\begin{array}{ll}\frac{1}{R}\sum _{\bar{n}}{\left[{\Psi }_{\bar{n}}^{(1)}({\mathbf{r}}^{\prime })\right]}^{\ast }{\Psi }_{\bar{n}}^{(3)}(\mathbf{r}),& \text{mbox}ifr>r^{\prime }\\ \frac{1}{R}\sum _{\bar{n}}{\left[{\Psi }_{\bar{n}}^{(3)}({\mathbf{r}}^{\prime })\right]}^{\ast }{\Psi }_{\bar{n}}^{(1)}(\mathbf{r}),& \text{mbox}ifr<r^{\prime }\end{array}$$

$${\Phi }_{\mathrm{ext}}(\mathbf{r})=$$

$${\Phi }_{\mathrm{in}}(\mathbf{r})=$$

$${\Phi }_{\mathrm{sca}}(\mathbf{r})=$$

$$\displaystyle\begin{array}[]{rl}-\tilde{q}_{\bar{n}}=&(\epsilon-1)\sum_{\bar{k}}~{}\tilde{L}^{11}_{\bar{n}\bar{k}}~{}\tilde{d}_{\bar{k}},\\[14.22636pt] \tilde{b}_{\bar{n}}-\tilde{d}_{\bar{n}}=&(\epsilon-1)\sum_{\bar{k}}~{}\tilde{L}^{31}_{\bar{n}\bar{k}}~{}\tilde{d}_{\bar{k}},\end{array}$$

$${\tilde{L}}_{\bar{n}\bar{k}}^{11}=\frac{1}{R}{\int }_S{\left[{\Psi }_{\bar{n}}^{(1)}({\mathbf{r}}^{\prime })\right]}^{\ast }\frac{\partial {\Psi }_{\bar{k}}^{(1)}}{\partial n^{\prime }}({\mathbf{r}}^{\prime })d{S}^{\prime },$$

$${\tilde{L}}_{\bar{n}\bar{k}}^{31}=\frac{1}{R}{\int }_S{\left[{\Psi }_{\bar{n}}^{(3)}({\mathbf{r}}^{\prime })\right]}^{\ast }\frac{\partial {\Psi }_{\bar{k}}^{(1)}}{\partial n^{\prime }}({\mathbf{r}}^{\prime })d{S}^{\prime }.$$

$$\tilde{\mathbf{q}}=-\tilde{\mathbf{P}}\tilde{\mathbf{d}},$$

$$\tilde{\mathbf{P}}$$

$$\tilde{\mathbf{Q}}$$

$$\tilde{\mathbf{T}}$$

$${\mathbf{N}}_{nm}^{(1)}$$

$$\text{vspace}.2cm{\mathbf{N}}_{nm}^{(3)}$$

$$B_n=\frac{1}{(2n-1)!!}\sqrt{\frac{(n+1)}{n(2n+1)}}.$$

$${\tilde{b}}_{\bar{n}}$$

$${\tilde{d}}_{\bar{n}}$$

$${\tilde{q}}_{\bar{n}}$$

$$T_{\bar{n}\bar{k}}^{22}(\lambda \to \infty )$$

$$P_{\bar{n}\bar{k}}^{22}(\lambda \to \infty )$$

$$Q_{\bar{n}\bar{k}}^{22}(\lambda \to \infty )$$

$${\tilde{T}}_{\bar{n}\bar{k}}\propto R^{-(n+k+1)},$$

$${\tilde{P}}_{\bar{n}\bar{k}}\propto R^{-(n+k+1)},$$

$${\tilde{Q}}_{\bar{n}\bar{k}}\propto R^{(n-k)}.$$

$$\mathbf{n}\mathrm{d}S=\left[r{\mathbf{e}}_r-r_{\theta }{\mathbf{e}}_{\theta }\right]r\sin \theta \mathrm{d}\theta \mathrm{d}\varphi$$

$$\frac{\partial {\Psi }_{km}^{(1)}}{\partial n}\mathrm{d}S$$

$$d_n^m(\theta )=$$

$${\tau }_n^m(\theta )=$$

Full text

1 Introduction

In electromagnetic scattering, it is often convenient to express the electric field in terms of multipoles, e.g.

[TABLE]

where $\mathbf{M}_{nm}$ and $\mathbf{N}_{nm}$ are the electric and magnetic multipolar fields, also called vector spherical wavefunctions, and form a complete basis of divergence-less solutions of Helmoltz equation found using separation of variables in spherical coordinates [1, 2, 3, 4]. The coefficients of the series expansion for the scattered field depend linearly on those of the incident field, and this relation defines the transition $T$-matrix [3]. Within the extended boundary condition method (EBCM) [5, 3], also called null-field method [6] or $T$-matrix method, the $T$-matrix is computed from the division of two matrices whose matrix elements are given by surface integrals on the particle surface. This approach was developed by Waterman in the 1960-70s [7, 8, 9] and remains one of the most powerful techniques for the study of electromagnetic scattering by particles. It has been extensively studied and applied in many areas to compute the scattering properties of particles of arbitrary shapes [10] with a particular emphasis on axisymmetric particles.

Only a few studies have however studied its potential use to solve electrostatic problems involving dielectric particles. One approach is to take the long-wavelength limit of the standard EBCM/$T$-matrix formalism. This was discussed by Waterman in the special two-dimensional case of infinite cylinders [11]. More recently, Farafonov and co-workers have developed an electrostatic (ES) equivalent of the $T$-matrix/EBCM formulation [12, 13], hereafter referred to as ES-EBCM. They then exploited it to derive asymptotic expansions (for large multipole orders) for the matrix elements arising from the EBCM formalism and to study the validity of the Rayleigh hypothesis for electromagnetic scattering by a particle. The hypothesis states that the scattered field expansion converges everywhere outside the particle [3, 4, 14]. While it is generally accepted that the Rayleigh hypothesis is not always valid, its range of validity and its link to the singularity of the solutions are still the subject of investigations [4, 15, 14, 13].

In this work we show that analytic expressions can in fact be obtained for all ES-EBCM matrix elements in the special case of spheroids. For this we use expansions of solid spherical harmonics in terms of solid spheroidal harmonics and vice-versa. The analytic $T$-matrix extends the previously-found asymptotic behavior of the matrix elements. It moreover provides an analytic illustration of the region of validity of the Rayleigh hypothesis and its link to the analytical continuation of the series solution.

The manuscript is organized as follows. In Sec. 2, we re-derive a fully self-contained ES-EBCM formulation, and establish the link to the EBCM quantitatively. This allows us to discuss the symmetry properties of the ES-EBCM matrices from their EBCM counterparts. We then focus on the case of spheroids, first using the ES-EBCM integrals (Secs. 3). A simple expression is derived for the first column of the $T$-matrix, which allows us to discuss the Rayleigh hypothesis and its connection with the analytic continuation of the series of the scattered field. In Sec. 4 we then use the separation of variable method in spheroidal coordinates, and relationships between spherical and spheroidal harmonics, to obtain analytic expressions for all matrix elements.

2 The EBCM for the electrostatics problem

As described in Refs. [12, 13], the EBCM for electrostatics (ES-EBCM) can be derived from scratch following the same approach as for electromagnetic scattering. For practical reasons, the basis functions are different to the long-wavelength limit of the spherical wave-functions, which makes quantitative connections between the two approaches cumbersome. Here we will therefore use different notations and conventions to previous work [12, 13] to simplify the link to the full wave EBCM. We summarize below the main definitions and derivations of the ES-EBCM and provide this quantitative link with the EBCM.

We consider a dielectric particle embedded in a medium with dielectric constants $\epsilon_{2}$ and $\epsilon_{1}$, respectively, and an external field $\mathbf{E}_{\mathrm{ext}}(\mathbf{r})$. The corresponding potential is $\Phi_{\mathrm{ext}}(\mathbf{r})$, and we denote $\Phi_{\mathrm{out}}(\mathbf{r})$ and $\Phi_{\mathrm{in}}(\mathbf{r})$ the (total) potentials outside and inside the particle.

2.1 Surface-integral formulation

The electrostatics problem can be formulated in terms of surface integral equations [16]. For two points $\mathbf{r}$ and $\mathbf{r^{\prime}}$, the free space Green’s function is

[TABLE]

Using Laplace’s equation, $\nabla^{2}\Phi_{\mathrm{in}}=0$, and Green’s second identity with $\Phi_{\mathrm{in}}$ and $G$ on the interior volume, we obtain:

[TABLE]

where $\partial/\partial n^{\prime}_{-}$ is the derivative just inside the surface, with respect to surface normal (and $\partial/\partial n^{\prime}_{+}$ is the derivative just outside).

For the exterior volume $\bar{V}$, we note that $\Phi_{\mathrm{out}}$ does not strictly satisfy Laplace’s equation everywhere because of sources $\rho_{\infty}$ at infinity associated with the external potential, i.e. $\nabla^{2}\Phi_{\mathrm{ext}}=-\rho_{\infty}/\epsilon_{0}$ and $\Phi_{\mathrm{ext}}=\int_{\bar{V}}G\rho_{\infty}/\epsilon_{0}dV^{\prime}$. Using Green’s second identity with $\Phi_{\mathrm{out}}$ and $G$ on the exterior volume, we obtain (note the change of sign to retain the normal pointing outward):

[TABLE]

where $\Phi_{\mathrm{out}}=\Phi_{\mathrm{ext}}+\Phi_{\mathrm{sca}}$, with $\Phi_{\mathrm{sca}}$ the “scattered” potential. 111The denomination “scattered” potential is aligned with the full-wave treatment, but the term “reflected” potential would be more suited to a purely electrostatics problem. This expression is analogous to the null-field equations of the EBCM.

Using the boundary conditions on the surface: $\Phi_{\mathrm{out}}=\Phi_{\mathrm{in}}$ and $\epsilon_{1}\partial\Phi_{\mathrm{out}}/\partial n_{+}=\epsilon_{2}\partial\Phi_{\mathrm{in}}/\partial n_{-}$, we can combine the above equations to eliminate the Green’s function derivative (note that $\partial G/\partial n^{\prime}_{+}=\partial G/\partial n^{\prime}_{-}$) and get:

[TABLE]

where $\epsilon=\epsilon_{2}/\epsilon_{1}$ is the relative dielectric constant.

2.2 Expansion in terms of solid spherical harmonics

The next step is to expand Green’s function and the potentials onto the basis of interest, here the standard solutions of Laplace’s equations in spherical coordinates, chosen here as:

[TABLE]

We include the arbitrary length $R$ so that the basis functions are dimensionless. Without it, the matrix elements would have dimensions dependent on the row and column, which is not ideal. The full-wave EBCM basis-functions are also commonly defined as dimensionless [3]. $R$ may be set as a characteristic length of the problem, although physical predictions are independent of the choice of basis functions and therefore independent of $R$.

Our choice of normalization is slightly different to that used by Farafonov et al., and will simplify the connection with the full-wave EBCM and with the approach using spheroidal coordinates in Sec. 4. These basis functions satisfy the following orthogonality relations on a sphere of radius $R$:

[TABLE]

and Green’s function can be expanded as [17]:

[TABLE]

where the combined index notation $\bar{n}$ is equivalent to both indices $n,m$ (with $|m|\leq n$).

The external and internal potentials can be written as series in terms of $\Psi^{(1)}_{nm}$, which are regular at the origin, while the scattered potential is expanded in terms of $\Psi^{(3)}_{nm}$, which vanish at infinity. We here use notation in line with the electromagnetic treatment as presented in Ref. [3], explicitly:

[TABLE]

2.3 Derivation of the ES-EBCM formula

To solve the problem (determine $\Phi_{\mathrm{sca}}$ and $\Phi_{\mathrm{in}}$), we substitute the series expansions of the potentials and of Green’s function into the integral equations (Eqs. 5) and use the orthogonality of the basis functions to obtain

[TABLE]

where (similarly to Ref. [12]):

[TABLE]

The matrices $\mathbf{L}^{11}$ and $\mathbf{L}^{31}$ obtained by Farafonov [12, 13] differ from these by some normalization and dimensional coefficients due to different basis functions.

Eqs. 16 form an infinite-dimensional linear system between the expansion coefficients, and just as with the EBCM, these can be written succinctly using matrix notations:

[TABLE]

where $\tilde{\mathbf{b}},\tilde{\mathbf{d}}$, $\tilde{\mathbf{q}}$ are vectors containing elements $\tilde{b}_{\bar{n}}$, $\tilde{d}_{\bar{n}}$, $\tilde{q}_{\bar{n}}$, and $\tilde{\mathbf{P}}$, $\tilde{\mathbf{Q}}$, $\tilde{\mathbf{T}}$ are infinite-dimensional matrices given by:

[TABLE]

Physically, $\tilde{\mathbf{P}}$ is the coupling matrix for the multipolar component of order $\bar{n}$ of the field created outside by a multipolar field of order $\bar{k}$ inside. $\tilde{\mathbf{Q}}^{-1}$ is the coupling matrix for the multipolar component of order $\bar{n}$ of the field created inside the particle when excited by an external multipolar field of order $\bar{k}$. It is interesting to note that the matrices $\tilde{\mathbf{L}}^{11}$ and $\tilde{\mathbf{L}}^{31}$ are independent of permitivity and depend on shape only. In fact, the matrix $\tilde{\mathbf{L}}^{31}$ is a generalization of the depolarization factor introduced in the standard electrostatics solution of the ellipsoid [18].

We also note that Green’s second identity on $[\Psi_{\bar{n}}^{(1)}]^{*}$ and $\Psi^{(1)}_{\bar{k}}$ implies that the matrix $\tilde{\mathbf{L}}^{11}$ is Hermitian, i.e. $\tilde{L}^{11}_{\bar{n}\bar{k}}=[\tilde{L}^{11}_{\bar{k}\bar{n}}]^{*}$.

2.4 Link with the standard EBCM formulation

The connection between the ES-EBCM and EBCM formulations was discussed in Ref. [19], but with slightly different definitions and only for the $T$-matrix, so for completeness we here adapt and expand the main results. We use the standard EBCM notations of Ref. [3] except for the $\mathbf{RgQ}$ matrix denoted $\mathbf{P}$ here. The link can be established by considering the long-wavelength limit of the EBCM. We do not here consider the magnetostatic limit or the limit of the coupling between magnetic and electric multipoles, so only the block ${\mathbf{T}}^{22}$ relating to electric/electric multipole coupling needs to be considered. In the long wavelength limit, the electric multipole fields reduce to [19]:222Note that there is an error for ${\mathbf{N}}^{(3)}_{nm}$ in Eq. (27) of Ref. [19] as the exponent for $k$ not should be $n+1$, but $n+2$ (as corrected here).

[TABLE]

where $\mathbf{N}_{nm}^{(1)}$ and $\mathbf{N}_{nm}^{(3)}$ are the regular and irregular vector spherical wavefunctions for electric multipoles, $k$ is the wavenumber ($k_{2}$ inside the particle, $k_{1}$ outside), and we have defined for convenience:

[TABLE]

Comparing the general expansions of the incident and scattered electric fields and potentials via $\mathbf{E}=-\nabla\Phi$, and using 21 and 22 we derive:

[TABLE]

where $s=k_{2}/k_{1}=sqrt{\epsilon}$ is the relative refractive index of the particle, and $b_{\bar{n}}$, $d_{\bar{n}}$, $q_{\bar{n}}$ are the expansion coefficients of the electric multipole component of the incident, internal and scattered electric fields in the EBCM [3].

From those, we deduce the relations between the matrices in the $\lambda\rightarrow\infty$ limit:

[TABLE]

Again the choice of $R$ does not affect the long-wavelength limit of the EBCM. This can be checked from Eqs. 33 and 20 which show that the matrix elements in the ES-EBCM depend on $R$ as follows:

[TABLE]

2.5 Axisymmetric particles

Similarly to the EBCM, the surface integrals $\mathbf{L}^{11}$ and $\mathbf{L}^{31}$ simplify substantially for particles with symmetry of rotation (around the $z$-axis). One major simplification is that only integrals with the same projected angular momentum ($m$) are non-zero, since the integrals over $\phi$ are zero unless $m=m^{\prime}$. $m$ can then be taken as a fixed parameter (and is sometimes omitted from the notation for simplicity). The surface integrals defining $\tilde{\mathbf{L}}^{11}$ and $\tilde{\mathbf{L}}^{31}$ reduce to integrals over $\theta$ and are all real, which implies that $\tilde{\mathbf{L}}^{11}$ is symmetric: $\tilde{L}^{11}_{\bar{n}\bar{k}}=\tilde{L}^{11}_{\bar{k}\bar{n}}$. The matrix elements in $\mathbf{P}$, $\mathbf{Q}$, $\mathbf{T}$, are then also real if $\epsilon$ is real. In the EM case, the $\mathbf{T}^{22}$ matrix is symmetric by virtue of optical reciprocity [3]. From the link between the ECBM and ES-EBCM (Eq. 27) we deduce that the matrix $\tilde{\mathbf{T}}$ is also symmetric (even though it is not obvious from Eq. 20).

The integral expression can be further simplified by describing the particle surface in spherical coordinates by an equation of the form $r(\theta)$. Denoting $r_{\theta}=\mathrm{d}r/\mathrm{d}\theta$, the surface normal and elementary area are

[TABLE]

from which we deduce (using $\frac{\partial}{\partial n}=\mathbf{n}\cdot\nabla$)

[TABLE]

In line with the notation of Refs. [3, 20], we write333This expression assumes that $P_{n}^{m}(\cos\theta)$ is defined with the Condon-Shortley phase $(-1)^{m}$.:

[TABLE]

where $(-)^{m}$ is shorthand for $(-1)^{m}$. We can then deduce after a few manipulations:

[TABLE]

Note that the symmetry of $\tilde{L}^{11}_{nk}$ can be made more obvious using integrations by part to obtain:

[TABLE]

Note also that all matrices can be computed for $m\geq 0$ only, and the values for $m<0$ derive from

[TABLE]

and identical relations for the others matrices.

3 Spheroidal particles in the ES-EBCM formalism

As highlighted in Refs. [20, 21, 22], spheroidal particles are a special case in EM scattering in the context of the $T$-matrix method and this is also the case in electrostatics. The electrostatics case is even more special as there exists a full analytical solution of the problem with the separation of variable method. This provides a means to find analytic expressions for the entire $\mathbf{T}$, $\mathbf{P}$, and $\mathbf{Q}$ matrices (see sec. 4). For clarity, we focus the discussion in this section on the simple case of a spheroidal dielectric particle in a uniform external field. The general case will be presented in the next section.

We consider a prolate spheroid of semi-axes $a$ along $x,y$ and $c$ along $z$ ($c>a$). It has half-focal length $f=\sqrt{c^{2}-a^{2}}$ and eccentricity $e=f/c$. The surface is then defined in spherical coordinates as:

[TABLE]

3.1 Special properties

Because of the reflection symmetry with respect to the $z=0$ plane ($\theta\rightarrow\pi-\theta$), half of the integrals are zero, namely:

[TABLE]

Moreover, it was shown in the EM case [20] that the dominant terms in some of the other integrals (with $n+k$ even) vanish for spheroids, causing serious numerical issues. In electrostatics, the situation is much simpler, as the corresponding integrals vanish completely [20, 12]. As a result, $\tilde{\mathbf{Q}}^{m}$ is upper triangular for spheroids:

[TABLE]

which can be deduced from the identities proved in Ref. [20]. One immediate consequence is that for a uniform external field (which has $n=1$ only in its series expansion), the inside field is also uniform (it also has $n=1$ only in its expansion since $\tilde{\mathbf{Q}}^{-1}$ must also be upper triangular). This is related to the Eshelby conjecture [23, 24], which states that ellipsoids (of which spheroids are a special case) are the only geometry for which the inside field is uniform when placed in an external uniform field, so we may also conjecture that ellipsoids are the only geometry for which $\tilde{\mathbf{Q}}$ is upper triangular.

3.2 Integral approach for a uniform field

Following the standard EBCM approach, we could calculate the matrix $\tilde{\mathbf{T}}$ from the surface integrals in Eqs. 33 and the expressions in Eqs. 20. This method is however cumbersome so we apply it only for illustration to the elements $\tilde{T}_{n1}^{m=0}$, which are sufficient to solve the problem for a uniform external field along the $z$-axis $\mathbf{E}=E_{0}\mathbf{\hat{z}}$. The expansion of the external potential has only one term: $\tilde{b}_{1}^{0}=-\sqrt{4\pi}RE_{0}$. Since $\tilde{\mathbf{Q}}$ is upper triangular, $\tilde{\mathbf{R}}=\tilde{\mathbf{Q}}^{-1}$ is also upper triangular and $\tilde{R}_{11}^{0}=\left(\tilde{Q}_{11}^{0}\right)^{-1}$, consequently

[TABLE]

and $\tilde{d}_{n}^{0}=0$ for $n\geq 2$ (the internal field is uniform, as previously discussed). Moreover,

[TABLE]

and the “scattered” field expansion coefficients are given by $\tilde{q}_{n}^{0}=\tilde{T}_{n1}^{0}\tilde{b}_{1}^{0}$.

We must now calculate the two $L$-integrals. Firstly:

[TABLE]

This is the standard depolarization factor along the long axis of a prolate spheroid, $L_{z}$ [18, 2], which and can be expressed explicitly as

[TABLE]

The second integrals, $\tilde{L}_{n1}^{11,m=0}$, can be simplified to:

[TABLE]

Because of the reflection symmetry, this integrates to zero for $n$ even. For $n$ odd, accurate numerical evaluation of these integrals is difficult because of oscillations and cancellations. However it is possible to find a very simple analytic solution:

[TABLE]

This result is derived in A using the integral expression and will also result from the derivation presented in the next section. We deduce that:

[TABLE]

The scattered potential expansion is then (note that it is $R$-independent as expected):

[TABLE]

where we have introduced the standard dipolar polarizability $\alpha_{zz}$ of the prolate spheroid along its long axis [18, 2, 25]:

[TABLE]

3.3 Revisiting the Rayleigh hypothesis

Thanks to this simple analytic expression (Eq. 46), it is straightforward to infer that the series diverges for $r<f$ and converges for $r>f$. This boundary marks the range of validity of the Rayleigh hypothesis for a prolate spheroid in the ES-EBCM framework. It is the same boundary as found using asymptotic expressions for the coefficients [12, 13], and as found in full-wave EBCM calculations [14]. This is expected given the link between the ES-EBCM and the long-wavelength limit of the EBCM. The advantage of our analytic expression is that it provides a vivid illustration of the mathematical concepts involved in discussing the Rayleigh hypothesis, in particular its link to the analytic continuation of the solution, on which we now elaborate.

If we consider a high aspect ratio spheroid with $c/a>\sqrt{2}$, there is a region close to the surface (away from the tips), where $r<f$ and the scattered field expansion diverges. Despite this, the ES-EBCM (or EBCM) remains valid as it does not rely on the Rayleigh hypothesis and the expansion in Eq. 46 remains correct where it converges – it only fails in some regions near the surface. The reason for this failure is that the spherical basis chosen to expand the scattered field is not well suited to the intrinsic singularity of the solution. To see this explicitly, one can find an analytic expression for the series expansion in Eq. 46 (derived in B):

[TABLE]

where $r^{\prime},z^{\prime}$ and $r^{\prime\prime},z^{\prime\prime}$ refer to the coordinates in the offset frames centered on the focal points F’ (at ${\mathbf{r}}=f\hat{\mathbf{z}}$) and F” (at ${\mathbf{r}}=-f\hat{\mathbf{z}}$), or explicitly:

[TABLE]

This is the analytic continuation of Eq. 46; it is singular only on the $z$ axis between the two foci of the spheroid. The basis of solid spherical harmonics was not well suited for such a solution because a series of these harmonics always has a spherical boundary of convergence, and such an expansion (Eq. 46) therefore increased the singular region from a line to the smallest sphere containing that line (the sphere $r=f$). Note nevertheless that since the analytic continuation of a series is unique, the series in Eq. 46 already contains, at least formally, all the information about the scattered field, even in the regions where it diverges. In general however, finding an analytic continuation is a difficult problem, and other methods must be used to calculate the scattered field near the surface [14]. These arguments have been laid out before from an abstract point of view [26, 27, 28, 29, 30, 31, 4], but we believe the simple analytic expressions derived here provide an insightful illustration.

4 Full analytic solution of the ES-EBCM formalism for prolate spheroids

In the case of the electrostatics problem for spheroids, the general solution may also be found using separation of variables in spheroidal coordinates. Those solutions can then be expressed in the basis of solid spherical harmonics to find analytic expressions for all the matrix elements of the ES-EBCM.

4.1 Separation of variables method

We use prolate spheroidal coordinates $(\xi,\eta,\phi)$ defined as follows [1]:

[TABLE]

where $r^{\prime}$, $r^{\prime\prime}$ were defined in Eq. 49 and $\phi$ is the same as in spherical coordinates. Surfaces of constant $\xi$ represent prolate spheroids with focal length $2f$. We will consider a dielectric prolate spheroid of semi-axes $a$ and $c$ ($c>a$) as in the previous section in a general external field. The half focal length is $f=\sqrt{c^{2}-a^{2}}$ and the eccentricity $e=f/c$. In spheroidal coordinates, the surface is at $\xi=\xi_{0}=1/e$.

For the most general external potential, this problem can be solved using solid spheroidal harmonic bases by expanding the potentials as:

[TABLE]

where $Q_{n}^{m}$ are the associated Legendre functions of the second kind.

The boundary conditions are:

[TABLE]

Applying these and using the fact that $P_{n}^{m}(\eta)$ form a basis, we have

[TABLE]

where $\partial P_{n}^{m}$ and $\partial Q_{n}^{m}$ denote the derivatives of $P_{n}^{m}$ and $Q_{n}^{m}$. From this, we derive

[TABLE]

where

[TABLE]

We have introduced the external $\Upsilon_{n}^{m}$ and internal $\Gamma_{n}^{m}$ spheroidal susceptibilities, which are similar to the multipolar polarizabilities for a sphere but depend on $m$, as expected for a non-spherically-symmetric object. $\Upsilon_{n}^{m}$ is analogous to the $\tilde{\mathbf{T}}$ matrix: for a spheroid, $\tilde{\mathbf{T}}^{m}$ is diagonal in the spheroidal basis with diagonal elements $\Upsilon_{n}^{m}$. Similar observations can be made for the other matrices: $\tilde{\mathbf{P}}$, $\tilde{\mathbf{Q}}$, $\tilde{\mathbf{L}}^{11}$, $\tilde{\mathbf{L}}^{31}$. For example $(\Gamma_{n}^{m})^{-1}$ are the diagonal elements of $\tilde{\mathbf{Q}}^{m}$ and $-\Upsilon_{n}^{m}/\Gamma_{n}^{m}$ the diagonal elements of $\tilde{\mathbf{P}}$ in the spheroidal basis.

The spheroidal susceptibilities can be simplified by identifying the Wronskian of the associated Legendre equation (see Sec. 8.1.8 in Ref. [32]):

[TABLE]

Then we have

[TABLE]

We can also rearrange the solution as:

[TABLE]

where

[TABLE]

are the diagonal elements of $\tilde{\mathbf{L}}^{11}$, $\tilde{\mathbf{L}}^{31}$ expressed in the spheroidal basis.

4.2 Analytic expressions for all matrix elements

Using the simple relationships between the incident, scattered and internal fields in the spheroidal basis, and the expansions between spherical and spheroidal harmonics, we can derive the equivalent relationships in the spherical basis and therefore analytic expressions for all elements of $\tilde{\mathbf{P}}$, $\tilde{\mathbf{Q}}$, and $\tilde{\mathbf{T}}$. The spherical-spheroidal harmonic expansions [33, 34] are summarized in C. These express the spherical basis functions, $r^{n}P_{n}^{m}(\cos\theta)$ or $r^{-n-1}P_{n}^{m}(\cos\theta)$, as a sum of spheroidal basis functions, $P_{k}^{m}(\xi)P_{k}^{m}(\eta)$ or $Q_{k}^{m}(\xi)P_{k}^{m}(\eta)$, and vice versa. We set $R=f$ in the following to simplify the resulting expressions.

Our goal is to express $\tilde{q}_{n}^{m}$ in terms of $\tilde{b}_{n}^{m}$. We first express $\alpha_{n}^{m}$ in terms of $\tilde{b}_{n}^{m}$, i.e. derive the spheroidal expansion of the external field from its spherical expansion. Inserting Eq. 99 into the external potential spherical expansion (Eq. 11), rearranging the order of summation, and identifying the coefficients $\alpha_{n}^{m}$ from Eq. 51, we obtain:

[TABLE]

We then express $\tilde{q}_{n}^{m}$ in terms of $\gamma_{n}^{m}$ by substituting Eq. 101 into Eq. 53 and identifying the coefficients $\tilde{q}_{n}^{m}$ from Eq. 13

[TABLE]

The relationship between $\tilde{q}_{n}^{m}$ and $\tilde{b}_{n}^{m}$ is found by combining Eqs. 66, 68 with the relation $\gamma_{s}^{m}=\alpha_{s}^{m}\Upsilon_{s}^{m}$ (Eq. 56) and swapping the order of summation:

[TABLE]

where $\min(n,k)$ denotes the minimum of $n$ and $k$. The matrix elements of $\mathbf{\tilde{T}}$ can be readily identified as:

[TABLE]

where we have defined

[TABLE]

Note that $\Upsilon_{s}^{m}$ depends on the spheroid shape (via $\xi_{0}$) and its optical properties (via $\epsilon$).

For $\tilde{\mathbf{L}}^{11}$, we use Eqs. 62 and the derivation is almost identical to that of $\tilde{\mathbf{T}}$:

[TABLE]

In this form, the symmetry of $\tilde{\mathbf{T}}$ and $\tilde{\mathbf{L}}^{11}$ is clear. For $\tilde{\mathbf{L}}^{31}$, we need the inverse relation of Eq. 66, which is obtained using the expansion Eq. 98. We then use Eq. 63 to get:

[TABLE]

As previously observed, the $L$-matrices only depend on shape, not on $\epsilon$. $\tilde{\mathbf{P}}$ and $\tilde{\mathbf{Q}}$ are easily derived from Eqs. 20.

4.3 Special cases

It is interesting to consider a few special cases. Firstly for $m=0$:

[TABLE]

The latter equality arises when setting $e=1/\xi_{0}$ in Eq. 42. We also have

[TABLE]

From those we deduce

[TABLE]

For the $T$-matrix, we have from Eq. 70:

[TABLE]

which is the same as Eq. 45 found earlier.

For $m=1$, we have

[TABLE]

The latter is the depolarization factor of the spheroid along the $x,y$ directions: $L_{x}=L_{y}=(1-L_{z})/2$. From those, we have

[TABLE]

and Eq. 70 simplifies to

[TABLE]

Finally, for a general $n$ and $m$, Eq. 73 also reduces to a single term when $k=n$:

[TABLE]

Those diagonal elements of $\tilde{L}^{31}$ can be viewed as the generalized depolarization factors of the spheroid for multipole orders $n,m$.

4.4 Comparison with EBCM results

To assess the applicability of our analytic results, we compare the analytic matrix elements for spheroids to those obtained numerically using an accurate implementation of the EBCM method for spheroids [22] for which codes are available [35]. The relations established in Eq. 27 are used to compare the two methods. For particles much smaller than the wavelength, the analytic quasi-static solution is expected to provide a good approximation to the exact result. This is illustrated in Fig. 1 for three dielectric prolate spheroids all of aspect ratio 10 and long semi-axis 10, 20 and 30 nm, excited with 600 nm wavelength light. The agreement between the two methods increases as the size decreases, as expected. Note however that the other blocks of the EBCM $T$-matrix, $\mathbf{T}^{11},~{}\mathbf{T}^{12}$, and $\mathbf{T}^{21}$ contain elements that are not negligible in the quasi-static limit. So for the ES-EBCM to be applicable as an approximation to quasi-static (not just electrostatic) problems, these other blocks should also be taken into account.

5 Conclusion

In summary, this work provides a detailed investigation of the electrostatics limit of Waterman’s $T$-matrix formalism, in particular in the special case of spheroidal particles. We focused specifically on prolate spheroids but the results could be adapted to oblate spheroids. We believe that the analytic expressions that we derived will provide an insight into the fundamental properties of the $T$-matrix formalism. This was illustrated in the context of the Rayleigh Hypothesis and its connection to the analytic continuation of the solutions. The ES-EBCM formalism could be further developed to encompass the magnetostatic limit (relating to magnetic multipoles) and the long-wavelength limit of the coupling matrices between electric and magnetic multipoles. The ES-EBCM formalism could moreover be applied to the solution of other electrostatic problems.

Appendix A Proof of Eq. 44

We reproduce here the form of the integral

[TABLE]

and consider $n$ odd only (the integral is zero for $n$ even). We first substitute $r(\theta)$ from Eq. 36 and make the change of variable $x=\cos\theta$:

[TABLE]

Using the binomial expansion, this is rewritten as a series

[TABLE]

where

[TABLE]

Note that $n$ is odd and this is a generalized binomial coefficient with half-integer entries.

We then expand the powers $x^{2p+1}$ on the basis of Legendre polynomials as

[TABLE]

This is as a special case of Eq. 99, evaluated at $\eta=1$, $n=2p+1$, $m=0$, where the coefficients are

[TABLE]

Using the orthogonality of the Legendre polynomials, the integral can be expressed as

[TABLE]

$\alpha_{p,n}\beta_{n,p}$ can be shown to simplify to

[TABLE]

Re-indexing the sum with $k=p-(n-1)/2$:

[TABLE]

where we recognize the binomial series expansion for $\frac{1}{(1-e^{2})^{n/2+1}}$, which simplifies using $e=f/c$ and $1-e^{2}=a^{2}/c^{2}$ to

[TABLE]

As shown in Sec. 4, this can be alternatively derived from the spheroidal coordinate solution of the problem. Given the remarkable simplicity of the final expression for $\tilde{L}_{n1}^{11,m=0}$, one might have expected that a simpler proof would be possible.

Appendix B Analytic continuation of the multipole expansion for the scattered potential

We start from the series expansion given in the main text (rewritten here up to a constant factor):

[TABLE]

We first recognize that this expression has some resemblance with the generating function for the Legendre polynomials, namely:

[TABLE]

We therefore define

[TABLE]

from which we infer

[TABLE]

The integral for $G$ can be calculated analytically:

[TABLE]

Note that this result can also be obtained from Eq. 100 with $n=0$. The simpler expression given in Eq. 48 can be obtained by considering the offset coordinate frames centred at both foci (the offset coordinates are defined explicitly in Eqs. 49).

Appendix C Relationships between spherical and spheroidal solid harmonics

Below are the four possible relations between the regular and irregular spherical solid harmonics, and the regular and irregular prolate spheroidal harmonics. The azimuthal dependence $e^{\pm im\phi}$ is omitted since it is the same on both sides. The spheroidal coordinates are defined in Eq. 50. By setting $R=f$ in the definitions of our spherical basis functions, those can easily be rewritten in terms of $\Psi_{nm}^{(1)}$ and $\Psi_{nm}^{(3)}$. Derivations can be found in [33, 34].

[TABLE]

Acknowledgments

We acknowledge the support of the Royal Society of New Zealand (RSNZ) through a Marsden Grant (ECLR) and a Rutherford Discovery Fellowship (BA).