# Connection Formulae between Ellipsoidal and Spherical Harmonics with Applications to Fluid Dynamics and Electromagnetic Scattering

Authors: Michael Doschoris; Panayiotis Vafeas
Journal: Advances in Mathematical Physics
 (2015)

---

## Abstract

The environment of the ellipsoidal system, significantly more complex than the spherical one, provides the necessary settings for tackling boundary value problems  in anisotropic space. However, the theory of Lamé functions and ellipsoidal harmonics affiliated with the ellipsoidal system is rather complicated. A turning point  would reside in the existence of expressions interlacing these two different systems. Still, there is no simple way, if at all, to bridge the gap. The present paper  addresses this issue. We provide explicit formulas of specific ellipsoidal harmonics expressed in terms of their counterparts in the classical spherical system. These expressions  are then put into practice in the framework of physical applications.


---

## Body

## 1. Introduction

The ellipsoidal coordinate system, by its very nature, is demanding concealing numerous difficulties. The main reason can be associated with the acquisition of solutions for miscellaneous operators. Even in the case of the Laplacian, deriving the corresponding eigensolutions is a nontrivial task. The French engineer and mathematician Gabriel Lamé in the mid nineteenth century, following an ingenious argument, separated variables for the Laplace operator arriving at the functions, which carry nowadays his name. Taking the product of Lamé functions leads to the ellipsoidal harmonics.But the complications regarding the particular system do not end here. In contrast to the theory of spherical harmonics, only ellipsoidal harmonics of low order have been computed in closed form [1, 2]. Why? First of all, a recursive technique in order to generate Lamé functions does not exist. Although we know that Lamé functions are connected by three-term recurrence relations [3], to the authors knowledge no procedure calculating the corresponding coefficients has been proposed so far. This particular impediment forces us to undergo an involved algorithm from which the Lamé functions are determined. This essentially two-step operation requires the computation of the roots of polynomial functions allowing nontrivial solutions for the initial linear homogeneous systems. We note that the previously indicated algorithm can be applied analytically only for Lamé functions up to the seventh degree [2]. Higher-order terms demand computational implementations [4] introducing numerical instability, which in the sequence is transferred to the calculation of the corresponding Lamé function.The aforementioned hurdles could in theory be avoided on the assumption that the functions of Lamé and corresponding ellipsoidal harmonics would be able to be expressed in terms of Legendre functions and spherical harmonics, respectively. Although, in principle, the possibility exists, no general formulae connecting these functions are available. The absence of such relations is justified bearing in mind that ellipsoidal harmonics are not reducible in a straightforward and unique way to the corresponding spherical harmonics. Nonetheless, distinct ellipsoidal harmonics can be represented in terms of a finite set of spherical harmonics of degree equal or less in view of the initial ellipsoidal harmonic.The present communication aims towards this direction. The section following is devoted to a brief introduction of the peculiarities of the ellipsoidal system, Lamé functions, and ellipsoidal harmonics (missing details can be found in [1]). We then continue elaborating explicit relations connecting ellipsoidal harmonics with the spherical ones. Finally, Section 3, introducing two real-world problems, showcases the efficiency of the derived formulas.


## 2. Mathematical Background

### 2.1. The Confocal Ellipsoidal Coordinate System

The ellipsoidal coordinates(
ρ
,
μ
,
ν
) of each point (
x
1
,
x
2
,
x
3
) are
(1)
x
1
=
ρ
μ
ν
h
2
h
3
,
x
2
=
ρ
2
-
h
3
2
μ
2
-
h
3
2
h
3
2
-
ν
2
h
1
h
3
,
x
3
=
ρ
2
-
h
2
2
h
2
2
-
μ
2
h
2
2
-
ν
2
h
1
h
2
,
provided that 0
≤
ν
2
≤
h
3
2
≤
μ
2
≤
h
2
2
≤
ρ
2
<
+
∞. The coordinate ρ, which determines a family of confocal ellipsoids, is comparable to the radial variable r in spherical coordinates. On the other hand, the coordinates μ and ν specify a family of confocal hyperboloids of one and two sheets, respectively, and correlate to the angular variables θ and ϕ.In the ellipsoidal coordinate system each direction is unique providing a particular perception of anisotropic space. Since any direction exhibits its own character, the ellipsoidal coordinate system requires the introduction of a reference ellipsoid establishing the variations in angular dependence, a direct analogy to the unit sphere. This reference ellipsoid defined by(2)
∑
j
=
1
3
‍
x
j
a
j
2
=
1
is not one of a kind but must incorporate the physical reality at hand. The constants
(3)
h
1
2
=
a
2
2
-
a
3
2
,
h
2
2
=
a
1
2
-
a
3
2
,
h
3
2
=
a
1
2
-
a
2
2
are the squares of the semifocal distances and a
κ
,
κ
=
1,2
,
3 with 0
<
a
3
<
a
2
<
a
1
<
+
∞ fixed parameters determining the reference semiaxes.A decisive aspect when concerned with boundary value problems in ellipsoidal coordinates resides in the spectral decomposition of the Laplacian. Separating variables leads to three identical ordinary differential equations, known as Lamé’s equation. The corresponding solutions are the so-called Lamé functions of the first kindE
n
m, where n designates the set of nonnegative integers and m
=
1,2
,
…
,
2
n
+
1, as well as the matching second kind functions F
n
m.Analytically, the first, second, and third degree Lamé functions of the first kind are presented below. The variableu represents either variable ρ
∈
[
h
2
,
+
∞
), μ
∈
[
h
3
,
h
2
], or ν
∈
[
0
,
h
3
]. Therein,
(4)
E
0
1
(
u
)
=
1
,
E
1
κ
(
u
)
=
u
2
-
a
1
2
+
a
κ
2
,
κ
=
1,2
,
3
,
E
2
1
(
u
)
=
u
2
+
Λ
-
a
1
2
,
E
2
2
(
u
)
=
u
2
+
Λ
′
-
a
1
2
,
E
2
κ
+
l
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
-
a
1
2
+
a
l
2
,
E
2
κ
+
l
u
=
1
a
1
2
+
a
κ
2
κ
,
l
=
1,2
,
3
,
κ
≠
l
,
E
3
2
κ
-
1
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
+
Λ
κ
-
a
1
2
,
E
2
κ
+
l
u
11111111
=
1
a
1
2
+
a
κ
2
κ
=
1,2
,
3
,
E
3
2
κ
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
+
Λ
κ
′
-
a
1
2
,
E
2
κ
+
l
u
11111
1
=
1
a
1
2
+
a
κ
2
κ
=
1,2
,
3
,
E
3
7
(
u
)
=
u
u
2
-
h
3
2
u
2
-
h
2
2
,
where the constants Λ and Λ
′, as well as Λ
κ and Λ
κ
′,  κ
=
1,2
,
3, are given as
(5)
Λ
=
a
1
2
-
h
3
2
+
h
2
2
3
+
h
1
4
+
h
2
2
h
3
2
3
,
(6)
Λ
′
=
a
1
2
-
h
3
2
+
h
2
2
3
-
h
1
4
+
h
2
2
h
3
2
3
,
(7)
Λ
κ
=
a
1
2
-
h
2
2
+
(
1
+
2
δ
κ
1
)
h
3
2
+
h
κ
2
5
+
h
1
4
+
h
2
2
h
3
2
+
3
h
κ
2
5
,
κ
=
1,2
,
3
,
(8)
Λ
κ
′
=
a
1
2
-
h
2
2
+
(
1
+
2
δ
κ
1
)
h
3
2
+
h
κ
2
5
-
h
1
4
+
h
2
2
h
3
2
+
3
h
κ
2
5
,
κ
=
1,2
,
3
,
where δ represents the Kronecker delta. The above constants satisfy the following relations:
(9)
∑
j
=
1
3
‍
1
Λ
-
a
j
2
=
0
,
∑
j
=
1
3
‍
1
+
2
δ
κ
j
Λ
κ
-
a
j
2
=
0
,
κ
=
1,2
,
3
,
respectively.The corresponding Lamé functions of the second kind are(10)
F
n
m
(
u
)
=
E
n
m
(
u
)
I
n
m
(
u
)
,
for every n
=
0,1
,
2
,
… and m
=
1,2
,
…
,
2
n
+
1, where
(11)
I
n
m
(
u
)
=
∫
u
0
u
‍
d
t
E
n
m
(
t
)
2
t
2
-
h
3
2
t
2
-
h
2
2
is an elliptic integral.In view of problems where the boundary consists of a triaxial confocal ellipsoid, the product of two Lamé functions belonging to the same class defines the surface ellipsoidal harmonicsS
n
m; that is,
(12)
S
n
m
(
μ
,
ν
)
=
E
n
m
(
μ
)
E
n
m
(
ν
)
,
whereas
(13)
E
n
m
(
ρ
,
μ
,
ν
)
=
E
n
m
(
ρ
)
S
n
m
(
μ
,
ν
)
designate the interior ellipsoidal harmonics.On the other hand, the exterior ellipsoidal harmonics are specified as(14)
F
n
m
(
ρ
,
μ
,
ν
)
=
(
2
n
+
1
)
E
n
m
(
ρ
,
μ
,
ν
)
I
n
m
(
ρ
)
,
where
(15)
I
n
m
(
ρ
)
=
∫
ρ
+
∞
‍
d
t
E
n
m
(
t
)
2
t
2
-
h
3
2
t
2
-
h
2
2
,
ρ
≥
h
2
.


### 2.2. Connecting Ellipsoidal and Spherical Harmonics

We already mentioned in the introduction the paucity of general formulas associating ellipsoidal harmonics with the same degree or less spherical harmonics. Another way to comprehend this is the following. As the triaxial ellipsoid deteriorates to a sphere, the ellipsoidal harmonics reduce to the so-called spheroconal harmonics which are a form of spherical harmonics. The spheroconal system, which incorporates the radial coordinater of the spherical system with the coordinates of the ellipsoidal system that specify orientation over any ellipsoidal surface (
μ
,
ν
), is established on the same parameters a
κ, κ
=
1,2
,
3, thus preserving its ellipsoidal characteristic. Nevertheless, although it seems that a general framework cannot be established, it is possible to represent distinct ellipsoidal harmonics with reference to finite terms of spherical harmonics (see Figure 1 for an illustration).Figure 1
In order to obtain representations of ellipsoidal harmonics in terms of spherical harmonics and vice versa, one has to go through Cartesian coordinates (solid lines). A direct connection appears not to be feasible (dashed line).Forasmuch as the spherical harmonics form a complete set, any continuous function can be expanded in a series ofY
n
m
(
r
^
). Therefore,
(16)
E
p
q
r
=
∑
n
=
0
∞
‍
∑
m
=
-
n
n
‍
A
n
m
r
n
Y
n
m
r
^
,
p
≥
0
,
q
=
1,2
,
…
,
2
p
+
1
and the coefficients A
n
m depend solely on the reference semiaxes a
1
,
a
2
,
a
3.Provided that on the unit sphereS
2
(17)
∮
S
2
‍
Y
n
m
(
r
^
)
Y
¯
n
′
m
′
(
r
^
)
d
Ω
(
r
^
)
=
δ
n
n
′
δ
m
m
′
,
the coefficients of (16) are computed as
(18)
A
n
m
=
1
r
n
∮
S
2
‍
E
p
q
(
r
)
Y
¯
n
m
(
r
^
)
d
Ω
(
r
^
)
.
Equations (16) and (18) provide the backbone of the presented analysis.Employing the Cartesian form of the ellipsoidal harmonics for degree up to three, namely,(19)
E
0
1
(
r
)
=
1
E
1
κ
(
r
)
=
h
1
h
2
h
3
h
κ
x
κ
,
κ
=
1,2
,
3
,
E
2
1
(
r
)
=
L
∑
j
=
1
3
‍
x
j
2
Λ
-
a
j
2
+
1
,
E
2
2
(
r
)
=
L
′
∑
j
=
1
3
‍
x
j
2
Λ
′
-
a
j
2
+
1
,
E
2
6
-
κ
(
r
)
=
h
1
h
2
h
3
h
κ
x
1
x
2
x
3
x
κ
,
κ
=
1,2
,
3
,
E
2
2
κ
-
1
(
r
)
=
h
1
h
2
h
3
L
κ
x
κ
h
κ
∑
j
=
1
3
‍
x
j
2
Λ
κ
-
a
j
2
+
1
,
κ
=
1,2
,
3
,
E
2
2
κ
(
r
)
=
h
1
h
2
h
3
L
κ
′
x
κ
h
κ
∑
j
=
1
3
‍
x
j
2
Λ
κ
′
-
a
j
2
+
1
,
κ
=
1,2
,
3
,
E
3
7
(
r
)
=
h
1
2
h
2
2
h
3
2
x
1
x
2
x
3
,
provided that r
=
(
x
,
y
,
z
) and
(20)
L
=
Λ
-
a
1
2
Λ
-
a
2
2
Λ
-
a
3
2
,
L
′
=
Λ
′
-
a
1
2
Λ
′
-
a
2
2
Λ
′
-
a
3
2
,
L
κ
=
Λ
κ
-
a
1
2
Λ
κ
-
a
2
2
Λ
κ
-
a
3
2
,
κ
=
1,2
,
3
,
L
κ
′
=
Λ
κ
′
-
a
1
2
Λ
κ
′
-
a
2
2
Λ
κ
′
-
a
3
2
,
κ
=
1,2
,
3
,
where the constants Λ
,
Λ
′
,
Λ
κ
,
Λ
κ
′
,
κ
=
1,2
,
3 are given by (5)–(8), respectively, the following representations hold:
(21)
E
0
1
(
r
)
=
4
π
r
0
Y
0
0
(
r
^
)
,
E
1
1
(
r
)
=
h
2
h
3
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
,
E
1
2
(
r
)
=
i
h
1
h
3
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
,
E
1
3
(
r
)
=
2
h
1
h
2
π
3
r
Y
1
0
(
r
^
)
,
E
2
1
r
=
L
+
4
π
5
L
Λ
-
a
3
2
r
2
Y
2
0
(
r
^
)
E
2
1
r
=
+
2
π
15
L
1
Λ
-
a
1
2
-
1
Λ
-
a
2
2
E
2
1
r
=
×
r
2
Y
2
-
2
(
r
^
)
+
Y
2
2
(
r
^
)
,
E
2
2
r
=
L
′
+
4
π
5
L
′
Λ
′
-
a
3
2
r
2
Y
2
0
(
r
^
)
E
2
1
r
=
+
2
π
15
L
′
1
Λ
′
-
a
1
2
-
1
Λ
′
-
a
2
2
E
2
1
r
=
×
r
2
Y
2
-
2
(
r
^
)
+
Y
2
2
(
r
^
)
,
E
2
3
(
r
)
=
i
h
1
h
2
h
3
2
2
π
15
r
2
Y
2
-
2
(
r
^
)
-
Y
2
2
(
r
^
)
,
E
2
4
(
r
)
=
h
1
h
2
2
h
3
2
π
15
r
2
Y
2
-
1
(
r
^
)
-
Y
2
1
(
r
^
)
,
E
2
5
(
r
)
=
i
h
1
2
h
2
h
3
2
π
15
r
2
Y
2
-
1
(
r
^
)
+
Y
2
1
(
r
^
)
,
E
3
1
r
=
h
2
h
3
L
1
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
+
h
2
h
3
L
1
Λ
1
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
-
Y
3
1
(
r
^
)
+
h
2
h
3
L
1
π
35
1
Λ
1
-
a
1
2
-
1
Λ
1
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
-
Y
3
3
(
r
^
)
,
E
3
2
r
=
h
2
h
3
L
1
′
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
+
h
2
h
3
L
1
′
Λ
1
′
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
-
Y
3
1
(
r
^
)
+
h
2
h
3
L
1
′
π
35
1
Λ
1
′
-
a
1
2
-
1
Λ
1
′
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
-
Y
3
3
(
r
^
)
,
E
3
3
r
=
i
h
1
h
3
L
2
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
+
i
h
1
h
3
L
2
Λ
2
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
+
Y
3
1
(
r
^
)
+
i
h
1
h
3
L
2
π
35
1
Λ
2
-
a
1
2
-
1
Λ
2
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
+
Y
3
3
(
r
^
)
,
E
3
4
r
=
i
h
1
h
3
L
2
′
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
+
i
h
1
h
3
L
2
′
Λ
2
′
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
+
Y
3
1
(
r
^
)
+
i
h
1
h
3
L
2
′
π
35
1
Λ
2
′
-
a
1
2
-
1
Λ
2
′
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
+
Y
3
3
(
r
^
)
,
E
3
5
r
=
2
h
1
h
2
L
3
π
3
r
Y
1
0
(
r
^
)
+
2
h
1
h
2
L
3
Λ
3
-
a
3
2
π
7
r
3
Y
3
0
(
r
^
)
+
h
1
h
2
L
3
2
π
105
1
Λ
3
-
a
1
2
-
1
Λ
3
-
a
2
2
×
r
3
Y
3
-
2
(
r
^
)
+
Y
3
2
(
r
^
)
,
E
3
6
r
=
2
h
1
h
2
L
3
′
π
3
r
Y
1
0
(
r
^
)
+
2
h
1
h
2
L
3
′
Λ
3
′
-
a
3
2
π
7
r
3
Y
3
0
(
r
^
)
+
h
1
h
2
L
3
′
2
π
105
1
Λ
3
′
-
a
1
2
-
1
Λ
3
′
-
a
2
2
×
r
3
Y
3
-
2
(
r
^
)
+
Y
3
2
(
r
^
)
,
E
3
7
(
r
)
=
i
h
1
2
h
2
2
h
3
2
2
π
105
r
3
Y
3
-
2
(
r
^
)
-
Y
3
2
(
r
^
)
.The above relations are not hard to prove, considering that the first term of the product present in the integrand of (18) displays the general form x
1
κ
x
2
l
x
3
s where κ
+
l
+
s
=
n, implying that only terms of the same or lower degree n survive. In addition, switching to spherical coordinates and bearing in mind that the point (
x
1
,
x
2
,
x
3
) resides interior of the ellipsoid (2), give the desired results.On the other hand, in order to evaluate the exterior ellipsoidal harmonics provided via (14) and (15) one needs only to express the quadratic terms E
n
m
(
t
)
2 as a function of Legendre polynomials. This is easily done furnishing
(22)
E
0
1
t
2
=
P
0
(
t
)
=
1
,
E
1
1
t
2
=
2
3
P
2
(
t
)
+
1
3
,
E
1
2
t
2
=
2
3
P
2
(
t
)
+
1
3
-
h
3
2
,
E
1
3
(
t
)
2
=
2
3
P
2
(
t
)
+
1
3
-
h
2
2
,
E
2
1
t
2
=
8
35
P
4
(
t
)
+
4
1
7
+
Λ
-
a
1
2
3
P
2
(
t
)
+
1
5
+
2
3
Λ
-
a
1
2
+
Λ
-
a
1
2
2
,
E
2
2
t
2
=
8
35
P
4
(
t
)
+
4
1
7
+
Λ
′
-
a
1
2
3
P
2
(
t
)
+
1
5
+
2
3
Λ
′
-
a
1
2
+
Λ
′
-
a
1
2
2
,
E
2
3
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
3
2
3
P
2
(
t
)
+
1
5
-
h
3
2
3
,
E
2
4
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
2
2
3
P
2
(
t
)
+
1
5
-
h
2
2
3
,
E
2
5
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
2
2
+
h
3
2
3
P
2
(
t
)
+
1
5
-
h
2
2
+
h
3
2
3
+
h
2
2
h
3
2
,
E
3
1
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
2
5
(
Λ
1
-
a
1
2
)
P
4
(
t
)
+
2
5
21
+
4
7
Λ
1
-
a
1
2
+
1
3
Λ
1
-
a
1
2
2
P
2
(
t
)
+
1
7
+
2
5
(
Λ
1
-
a
1
2
)
+
1
3
Λ
1
-
a
1
2
2
,
E
3
2
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
2
5
(
Λ
1
′
-
a
1
2
)
P
4
(
t
)
+
2
5
21
+
4
7
(
Λ
1
′
-
a
1
2
)
+
1
3
Λ
1
′
-
a
1
2
2
P
2
(
t
)
+
1
7
+
2
5
(
Λ
1
′
-
a
1
2
)
+
1
3
Λ
1
′
-
a
1
2
2
,
E
3
3
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
2
-
a
1
2
)
-
h
3
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
2
-
a
1
2
)
-
h
3
2
+
1
3
Λ
2
-
a
1
2
Λ
2
-
a
1
2
-
2
h
3
2
P
2
t
+
1
7
+
1
5
2
Λ
2
-
a
1
2
-
h
3
2
+
1
3
Λ
2
-
a
1
2
1
7
×
(
Λ
2
-
a
1
2
)
-
2
h
3
2
-
h
3
2
Λ
2
-
a
1
2
2
,
E
3
4
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
+
1
3
(
Λ
2
′
-
a
1
2
)
(
Λ
2
′
-
a
1
2
)
-
2
h
3
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
+
1
3
(
Λ
2
′
-
a
1
2
)
1
3
×
(
Λ
2
′
-
a
1
2
)
-
2
h
3
2
-
h
3
2
Λ
2
′
-
a
1
2
2
,
E
3
5
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
3
-
a
1
2
)
-
h
2
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
3
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
-
a
1
2
)
(
Λ
3
-
a
1
2
)
-
2
h
2
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
3
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
-
a
1
2
)
1
3
×
(
Λ
3
-
a
1
2
)
-
2
h
2
2
-
h
2
2
Λ
3
-
a
1
2
2
,
E
3
6
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
′
-
a
1
2
)
(
Λ
3
′
-
a
1
2
)
-
2
h
2
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
′
-
a
1
2
)
1
7
×
(
Λ
3
′
-
a
1
2
)
-
2
h
2
2
-
h
2
2
Λ
3
′
-
a
1
2
2
,
E
3
7
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
-
1
5
(
h
2
2
+
h
3
2
)
P
4
(
t
)
+
2
5
21
-
2
7
(
h
2
2
+
h
3
2
)
+
1
3
(
h
2
2
h
3
2
)
P
2
(
t
)
+
1
7
-
1
5
(
h
2
2
+
h
3
2
)
+
1
3
(
h
2
2
h
3
2
)
.


## 2.1. The Confocal Ellipsoidal Coordinate System

The ellipsoidal coordinates(
ρ
,
μ
,
ν
) of each point (
x
1
,
x
2
,
x
3
) are
(1)
x
1
=
ρ
μ
ν
h
2
h
3
,
x
2
=
ρ
2
-
h
3
2
μ
2
-
h
3
2
h
3
2
-
ν
2
h
1
h
3
,
x
3
=
ρ
2
-
h
2
2
h
2
2
-
μ
2
h
2
2
-
ν
2
h
1
h
2
,
provided that 0
≤
ν
2
≤
h
3
2
≤
μ
2
≤
h
2
2
≤
ρ
2
<
+
∞. The coordinate ρ, which determines a family of confocal ellipsoids, is comparable to the radial variable r in spherical coordinates. On the other hand, the coordinates μ and ν specify a family of confocal hyperboloids of one and two sheets, respectively, and correlate to the angular variables θ and ϕ.In the ellipsoidal coordinate system each direction is unique providing a particular perception of anisotropic space. Since any direction exhibits its own character, the ellipsoidal coordinate system requires the introduction of a reference ellipsoid establishing the variations in angular dependence, a direct analogy to the unit sphere. This reference ellipsoid defined by(2)
∑
j
=
1
3
‍
x
j
a
j
2
=
1
is not one of a kind but must incorporate the physical reality at hand. The constants
(3)
h
1
2
=
a
2
2
-
a
3
2
,
h
2
2
=
a
1
2
-
a
3
2
,
h
3
2
=
a
1
2
-
a
2
2
are the squares of the semifocal distances and a
κ
,
κ
=
1,2
,
3 with 0
<
a
3
<
a
2
<
a
1
<
+
∞ fixed parameters determining the reference semiaxes.A decisive aspect when concerned with boundary value problems in ellipsoidal coordinates resides in the spectral decomposition of the Laplacian. Separating variables leads to three identical ordinary differential equations, known as Lamé’s equation. The corresponding solutions are the so-called Lamé functions of the first kindE
n
m, where n designates the set of nonnegative integers and m
=
1,2
,
…
,
2
n
+
1, as well as the matching second kind functions F
n
m.Analytically, the first, second, and third degree Lamé functions of the first kind are presented below. The variableu represents either variable ρ
∈
[
h
2
,
+
∞
), μ
∈
[
h
3
,
h
2
], or ν
∈
[
0
,
h
3
]. Therein,
(4)
E
0
1
(
u
)
=
1
,
E
1
κ
(
u
)
=
u
2
-
a
1
2
+
a
κ
2
,
κ
=
1,2
,
3
,
E
2
1
(
u
)
=
u
2
+
Λ
-
a
1
2
,
E
2
2
(
u
)
=
u
2
+
Λ
′
-
a
1
2
,
E
2
κ
+
l
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
-
a
1
2
+
a
l
2
,
E
2
κ
+
l
u
=
1
a
1
2
+
a
κ
2
κ
,
l
=
1,2
,
3
,
κ
≠
l
,
E
3
2
κ
-
1
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
+
Λ
κ
-
a
1
2
,
E
2
κ
+
l
u
11111111
=
1
a
1
2
+
a
κ
2
κ
=
1,2
,
3
,
E
3
2
κ
u
=
u
2
-
a
1
2
+
a
κ
2
u
2
+
Λ
κ
′
-
a
1
2
,
E
2
κ
+
l
u
11111
1
=
1
a
1
2
+
a
κ
2
κ
=
1,2
,
3
,
E
3
7
(
u
)
=
u
u
2
-
h
3
2
u
2
-
h
2
2
,
where the constants Λ and Λ
′, as well as Λ
κ and Λ
κ
′,  κ
=
1,2
,
3, are given as
(5)
Λ
=
a
1
2
-
h
3
2
+
h
2
2
3
+
h
1
4
+
h
2
2
h
3
2
3
,
(6)
Λ
′
=
a
1
2
-
h
3
2
+
h
2
2
3
-
h
1
4
+
h
2
2
h
3
2
3
,
(7)
Λ
κ
=
a
1
2
-
h
2
2
+
(
1
+
2
δ
κ
1
)
h
3
2
+
h
κ
2
5
+
h
1
4
+
h
2
2
h
3
2
+
3
h
κ
2
5
,
κ
=
1,2
,
3
,
(8)
Λ
κ
′
=
a
1
2
-
h
2
2
+
(
1
+
2
δ
κ
1
)
h
3
2
+
h
κ
2
5
-
h
1
4
+
h
2
2
h
3
2
+
3
h
κ
2
5
,
κ
=
1,2
,
3
,
where δ represents the Kronecker delta. The above constants satisfy the following relations:
(9)
∑
j
=
1
3
‍
1
Λ
-
a
j
2
=
0
,
∑
j
=
1
3
‍
1
+
2
δ
κ
j
Λ
κ
-
a
j
2
=
0
,
κ
=
1,2
,
3
,
respectively.The corresponding Lamé functions of the second kind are(10)
F
n
m
(
u
)
=
E
n
m
(
u
)
I
n
m
(
u
)
,
for every n
=
0,1
,
2
,
… and m
=
1,2
,
…
,
2
n
+
1, where
(11)
I
n
m
(
u
)
=
∫
u
0
u
‍
d
t
E
n
m
(
t
)
2
t
2
-
h
3
2
t
2
-
h
2
2
is an elliptic integral.In view of problems where the boundary consists of a triaxial confocal ellipsoid, the product of two Lamé functions belonging to the same class defines the surface ellipsoidal harmonicsS
n
m; that is,
(12)
S
n
m
(
μ
,
ν
)
=
E
n
m
(
μ
)
E
n
m
(
ν
)
,
whereas
(13)
E
n
m
(
ρ
,
μ
,
ν
)
=
E
n
m
(
ρ
)
S
n
m
(
μ
,
ν
)
designate the interior ellipsoidal harmonics.On the other hand, the exterior ellipsoidal harmonics are specified as(14)
F
n
m
(
ρ
,
μ
,
ν
)
=
(
2
n
+
1
)
E
n
m
(
ρ
,
μ
,
ν
)
I
n
m
(
ρ
)
,
where
(15)
I
n
m
(
ρ
)
=
∫
ρ
+
∞
‍
d
t
E
n
m
(
t
)
2
t
2
-
h
3
2
t
2
-
h
2
2
,
ρ
≥
h
2
.


## 2.2. Connecting Ellipsoidal and Spherical Harmonics

We already mentioned in the introduction the paucity of general formulas associating ellipsoidal harmonics with the same degree or less spherical harmonics. Another way to comprehend this is the following. As the triaxial ellipsoid deteriorates to a sphere, the ellipsoidal harmonics reduce to the so-called spheroconal harmonics which are a form of spherical harmonics. The spheroconal system, which incorporates the radial coordinater of the spherical system with the coordinates of the ellipsoidal system that specify orientation over any ellipsoidal surface (
μ
,
ν
), is established on the same parameters a
κ, κ
=
1,2
,
3, thus preserving its ellipsoidal characteristic. Nevertheless, although it seems that a general framework cannot be established, it is possible to represent distinct ellipsoidal harmonics with reference to finite terms of spherical harmonics (see Figure 1 for an illustration).Figure 1
In order to obtain representations of ellipsoidal harmonics in terms of spherical harmonics and vice versa, one has to go through Cartesian coordinates (solid lines). A direct connection appears not to be feasible (dashed line).Forasmuch as the spherical harmonics form a complete set, any continuous function can be expanded in a series ofY
n
m
(
r
^
). Therefore,
(16)
E
p
q
r
=
∑
n
=
0
∞
‍
∑
m
=
-
n
n
‍
A
n
m
r
n
Y
n
m
r
^
,
p
≥
0
,
q
=
1,2
,
…
,
2
p
+
1
and the coefficients A
n
m depend solely on the reference semiaxes a
1
,
a
2
,
a
3.Provided that on the unit sphereS
2
(17)
∮
S
2
‍
Y
n
m
(
r
^
)
Y
¯
n
′
m
′
(
r
^
)
d
Ω
(
r
^
)
=
δ
n
n
′
δ
m
m
′
,
the coefficients of (16) are computed as
(18)
A
n
m
=
1
r
n
∮
S
2
‍
E
p
q
(
r
)
Y
¯
n
m
(
r
^
)
d
Ω
(
r
^
)
.
Equations (16) and (18) provide the backbone of the presented analysis.Employing the Cartesian form of the ellipsoidal harmonics for degree up to three, namely,(19)
E
0
1
(
r
)
=
1
E
1
κ
(
r
)
=
h
1
h
2
h
3
h
κ
x
κ
,
κ
=
1,2
,
3
,
E
2
1
(
r
)
=
L
∑
j
=
1
3
‍
x
j
2
Λ
-
a
j
2
+
1
,
E
2
2
(
r
)
=
L
′
∑
j
=
1
3
‍
x
j
2
Λ
′
-
a
j
2
+
1
,
E
2
6
-
κ
(
r
)
=
h
1
h
2
h
3
h
κ
x
1
x
2
x
3
x
κ
,
κ
=
1,2
,
3
,
E
2
2
κ
-
1
(
r
)
=
h
1
h
2
h
3
L
κ
x
κ
h
κ
∑
j
=
1
3
‍
x
j
2
Λ
κ
-
a
j
2
+
1
,
κ
=
1,2
,
3
,
E
2
2
κ
(
r
)
=
h
1
h
2
h
3
L
κ
′
x
κ
h
κ
∑
j
=
1
3
‍
x
j
2
Λ
κ
′
-
a
j
2
+
1
,
κ
=
1,2
,
3
,
E
3
7
(
r
)
=
h
1
2
h
2
2
h
3
2
x
1
x
2
x
3
,
provided that r
=
(
x
,
y
,
z
) and
(20)
L
=
Λ
-
a
1
2
Λ
-
a
2
2
Λ
-
a
3
2
,
L
′
=
Λ
′
-
a
1
2
Λ
′
-
a
2
2
Λ
′
-
a
3
2
,
L
κ
=
Λ
κ
-
a
1
2
Λ
κ
-
a
2
2
Λ
κ
-
a
3
2
,
κ
=
1,2
,
3
,
L
κ
′
=
Λ
κ
′
-
a
1
2
Λ
κ
′
-
a
2
2
Λ
κ
′
-
a
3
2
,
κ
=
1,2
,
3
,
where the constants Λ
,
Λ
′
,
Λ
κ
,
Λ
κ
′
,
κ
=
1,2
,
3 are given by (5)–(8), respectively, the following representations hold:
(21)
E
0
1
(
r
)
=
4
π
r
0
Y
0
0
(
r
^
)
,
E
1
1
(
r
)
=
h
2
h
3
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
,
E
1
2
(
r
)
=
i
h
1
h
3
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
,
E
1
3
(
r
)
=
2
h
1
h
2
π
3
r
Y
1
0
(
r
^
)
,
E
2
1
r
=
L
+
4
π
5
L
Λ
-
a
3
2
r
2
Y
2
0
(
r
^
)
E
2
1
r
=
+
2
π
15
L
1
Λ
-
a
1
2
-
1
Λ
-
a
2
2
E
2
1
r
=
×
r
2
Y
2
-
2
(
r
^
)
+
Y
2
2
(
r
^
)
,
E
2
2
r
=
L
′
+
4
π
5
L
′
Λ
′
-
a
3
2
r
2
Y
2
0
(
r
^
)
E
2
1
r
=
+
2
π
15
L
′
1
Λ
′
-
a
1
2
-
1
Λ
′
-
a
2
2
E
2
1
r
=
×
r
2
Y
2
-
2
(
r
^
)
+
Y
2
2
(
r
^
)
,
E
2
3
(
r
)
=
i
h
1
h
2
h
3
2
2
π
15
r
2
Y
2
-
2
(
r
^
)
-
Y
2
2
(
r
^
)
,
E
2
4
(
r
)
=
h
1
h
2
2
h
3
2
π
15
r
2
Y
2
-
1
(
r
^
)
-
Y
2
1
(
r
^
)
,
E
2
5
(
r
)
=
i
h
1
2
h
2
h
3
2
π
15
r
2
Y
2
-
1
(
r
^
)
+
Y
2
1
(
r
^
)
,
E
3
1
r
=
h
2
h
3
L
1
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
+
h
2
h
3
L
1
Λ
1
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
-
Y
3
1
(
r
^
)
+
h
2
h
3
L
1
π
35
1
Λ
1
-
a
1
2
-
1
Λ
1
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
-
Y
3
3
(
r
^
)
,
E
3
2
r
=
h
2
h
3
L
1
′
2
π
3
r
Y
1
-
1
(
r
^
)
-
Y
1
1
(
r
^
)
+
h
2
h
3
L
1
′
Λ
1
′
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
-
Y
3
1
(
r
^
)
+
h
2
h
3
L
1
′
π
35
1
Λ
1
′
-
a
1
2
-
1
Λ
1
′
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
-
Y
3
3
(
r
^
)
,
E
3
3
r
=
i
h
1
h
3
L
2
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
+
i
h
1
h
3
L
2
Λ
2
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
+
Y
3
1
(
r
^
)
+
i
h
1
h
3
L
2
π
35
1
Λ
2
-
a
1
2
-
1
Λ
2
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
+
Y
3
3
(
r
^
)
,
E
3
4
r
=
i
h
1
h
3
L
2
′
2
π
3
r
Y
1
-
1
(
r
^
)
+
Y
1
1
(
r
^
)
+
i
h
1
h
3
L
2
′
Λ
2
′
-
a
3
2
π
21
r
3
Y
3
-
1
(
r
^
)
+
Y
3
1
(
r
^
)
+
i
h
1
h
3
L
2
′
π
35
1
Λ
2
′
-
a
1
2
-
1
Λ
2
′
-
a
2
2
×
r
3
Y
3
-
3
(
r
^
)
+
Y
3
3
(
r
^
)
,
E
3
5
r
=
2
h
1
h
2
L
3
π
3
r
Y
1
0
(
r
^
)
+
2
h
1
h
2
L
3
Λ
3
-
a
3
2
π
7
r
3
Y
3
0
(
r
^
)
+
h
1
h
2
L
3
2
π
105
1
Λ
3
-
a
1
2
-
1
Λ
3
-
a
2
2
×
r
3
Y
3
-
2
(
r
^
)
+
Y
3
2
(
r
^
)
,
E
3
6
r
=
2
h
1
h
2
L
3
′
π
3
r
Y
1
0
(
r
^
)
+
2
h
1
h
2
L
3
′
Λ
3
′
-
a
3
2
π
7
r
3
Y
3
0
(
r
^
)
+
h
1
h
2
L
3
′
2
π
105
1
Λ
3
′
-
a
1
2
-
1
Λ
3
′
-
a
2
2
×
r
3
Y
3
-
2
(
r
^
)
+
Y
3
2
(
r
^
)
,
E
3
7
(
r
)
=
i
h
1
2
h
2
2
h
3
2
2
π
105
r
3
Y
3
-
2
(
r
^
)
-
Y
3
2
(
r
^
)
.The above relations are not hard to prove, considering that the first term of the product present in the integrand of (18) displays the general form x
1
κ
x
2
l
x
3
s where κ
+
l
+
s
=
n, implying that only terms of the same or lower degree n survive. In addition, switching to spherical coordinates and bearing in mind that the point (
x
1
,
x
2
,
x
3
) resides interior of the ellipsoid (2), give the desired results.On the other hand, in order to evaluate the exterior ellipsoidal harmonics provided via (14) and (15) one needs only to express the quadratic terms E
n
m
(
t
)
2 as a function of Legendre polynomials. This is easily done furnishing
(22)
E
0
1
t
2
=
P
0
(
t
)
=
1
,
E
1
1
t
2
=
2
3
P
2
(
t
)
+
1
3
,
E
1
2
t
2
=
2
3
P
2
(
t
)
+
1
3
-
h
3
2
,
E
1
3
(
t
)
2
=
2
3
P
2
(
t
)
+
1
3
-
h
2
2
,
E
2
1
t
2
=
8
35
P
4
(
t
)
+
4
1
7
+
Λ
-
a
1
2
3
P
2
(
t
)
+
1
5
+
2
3
Λ
-
a
1
2
+
Λ
-
a
1
2
2
,
E
2
2
t
2
=
8
35
P
4
(
t
)
+
4
1
7
+
Λ
′
-
a
1
2
3
P
2
(
t
)
+
1
5
+
2
3
Λ
′
-
a
1
2
+
Λ
′
-
a
1
2
2
,
E
2
3
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
3
2
3
P
2
(
t
)
+
1
5
-
h
3
2
3
,
E
2
4
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
2
2
3
P
2
(
t
)
+
1
5
-
h
2
2
3
,
E
2
5
(
t
)
2
=
8
35
P
4
(
t
)
+
2
2
7
-
h
2
2
+
h
3
2
3
P
2
(
t
)
+
1
5
-
h
2
2
+
h
3
2
3
+
h
2
2
h
3
2
,
E
3
1
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
2
5
(
Λ
1
-
a
1
2
)
P
4
(
t
)
+
2
5
21
+
4
7
Λ
1
-
a
1
2
+
1
3
Λ
1
-
a
1
2
2
P
2
(
t
)
+
1
7
+
2
5
(
Λ
1
-
a
1
2
)
+
1
3
Λ
1
-
a
1
2
2
,
E
3
2
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
2
5
(
Λ
1
′
-
a
1
2
)
P
4
(
t
)
+
2
5
21
+
4
7
(
Λ
1
′
-
a
1
2
)
+
1
3
Λ
1
′
-
a
1
2
2
P
2
(
t
)
+
1
7
+
2
5
(
Λ
1
′
-
a
1
2
)
+
1
3
Λ
1
′
-
a
1
2
2
,
E
3
3
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
2
-
a
1
2
)
-
h
3
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
2
-
a
1
2
)
-
h
3
2
+
1
3
Λ
2
-
a
1
2
Λ
2
-
a
1
2
-
2
h
3
2
P
2
t
+
1
7
+
1
5
2
Λ
2
-
a
1
2
-
h
3
2
+
1
3
Λ
2
-
a
1
2
1
7
×
(
Λ
2
-
a
1
2
)
-
2
h
3
2
-
h
3
2
Λ
2
-
a
1
2
2
,
E
3
4
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
+
1
3
(
Λ
2
′
-
a
1
2
)
(
Λ
2
′
-
a
1
2
)
-
2
h
3
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
2
′
-
a
1
2
)
-
h
3
2
+
1
3
(
Λ
2
′
-
a
1
2
)
1
3
×
(
Λ
2
′
-
a
1
2
)
-
2
h
3
2
-
h
3
2
Λ
2
′
-
a
1
2
2
,
E
3
5
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
3
-
a
1
2
)
-
h
2
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
3
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
-
a
1
2
)
(
Λ
3
-
a
1
2
)
-
2
h
2
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
3
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
-
a
1
2
)
1
3
×
(
Λ
3
-
a
1
2
)
-
2
h
2
2
-
h
2
2
Λ
3
-
a
1
2
2
,
E
3
6
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
+
1
5
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
P
4
(
t
)
+
2
5
21
+
2
7
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
′
-
a
1
2
)
(
Λ
3
′
-
a
1
2
)
-
2
h
2
2
P
2
(
t
)
+
1
7
+
1
5
2
(
Λ
3
′
-
a
1
2
)
-
h
2
2
+
1
3
(
Λ
3
′
-
a
1
2
)
1
7
×
(
Λ
3
′
-
a
1
2
)
-
2
h
2
2
-
h
2
2
Λ
3
′
-
a
1
2
2
,
E
3
7
(
t
)
2
=
16
231
P
6
(
t
)
+
8
7
3
11
-
1
5
(
h
2
2
+
h
3
2
)
P
4
(
t
)
+
2
5
21
-
2
7
(
h
2
2
+
h
3
2
)
+
1
3
(
h
2
2
h
3
2
)
P
2
(
t
)
+
1
7
-
1
5
(
h
2
2
+
h
3
2
)
+
1
3
(
h
2
2
h
3
2
)
.


## 3. Applications

In order to demonstrate the practicality of the preceding relations two physical applications are illustrated. The first problem comprises the field of low-Reynolds number hydrodynamics, while the second one emerges from the area of low-frequency electromagnetic scattering.


### 3.1. Particle-in-Cell Models for Stokes Flow in Ellipsoidal Geometry

A dimensionless criterion, which determines the relative importance of inertial and viscous effects in fluid dynamics, is the Reynolds number (Re). The pseudosteady and nonaxisymmetric creeping flow (Re
≪
1) of an incompressible, viscous fluid with dynamic viscosity μ
f and mass density ρ
f is described by the well-known Stokes equations.Let us consider any smooth and bounded or unbounded three-dimensional domainV
R
3, depending on the physical problem, in terms of the position vector r
=
x
1
,
x
2
,
x
3. Therein, Stokes flow governing equations connect the biharmonic vector velocity v with the harmonic scalar total pressure field P via
(23)
μ
f
Δ
v
=
∇
P
,
∇
·
v
=
0
for every r
∈
V
R
3.The left-hand side of (23) states that, in creeping flow, the viscous force compensates the force caused by the pressure gradient, while the right-hand side of it secures the incompressibility of the fluid. The fluid total pressure is related to the thermodynamic pressure p through
(24)
P
=
p
+
ρ
f
g
h
,
where the contribution of the term ρ
f
g
h, with g being the acceleration of gravity, refers to the gravitational pressure force, corresponding to a height of reference h. Once the velocity field is obtained, the harmonic vorticity field ω is defined as
(25)
ω
=
1
2
∇
×
v
,
while the stress dyadic Π
~ is assumed to be
(26)
Π
~
=
-
p
I
~
+
μ
f
∇
⊗
v
+
∇
⊗
v
⊺
in terms of the unit dyadic I
~
=
∑
j
=
1
3
‍
x
^
j
⊗
x
^
j, where the symbols ⊗ and ⊺ denote juxtaposition and transposition, respectively.One of the important areas of applications concerns the construction of particle-in-cell models, which are useful in the development of simple but reliable analytical expressions for heat and mass transfer in swarms of particles in the case of concentrated suspensions. In applied type analysis it is not usually necessary to have detailed solution of the flow field over the entire swarm of particles taking into account the exact positions of the particles, since such solutions are cumbersome to use. Thus, the technique of cell models is adopted where the mathematical treatment of each problem is based on the assumption that a three-dimensional assemblage may be considered to consist of a number of identical unit cells. Each of these cells contains a particle surrounded by a fluid envelope, containing a volume of fluid sufficient to make the fractional void volume in the cell identical to that in the entire assemblage.With the aim of mathematical modeling of the particular physical problem, we inherit, chosen among others cited in original paper [5] from where the main results were drawn, the Happel cell model [6] in which both the particle and the outer envelope enjoy spherical symmetry. In view of the Happel-type model [5, 6], we consider a fluid-particle system consisting of any finite number of rigid particles of arbitrary shape. Introducing the particle-in-cell model, we examine the Stokes flow of one of the assemblages of particles neglecting the interaction with other particles or with the bounded walls of a container. Let S
a denote the surface of the particle of the swarm, which is solid, is moving with a known constant translational velocity U
=
(
U
1
,
U
2
,
U
3
) in an arbitrary direction, and is rotating, also arbitrarily, with a defined constant angular velocity Ω
=
(
Ω
1
,
Ω
2
,
Ω
3
). It lives within an otherwise quiescent fluid layer, which is confined by the outer surface denoted by S
b. Following the formulation of Happel [6], extended to three-dimensional flows, the necessary nonslip flow conditions
(27)
v
=
U
+
Ω
×
r
,
r
∈
S
a
are imposed on the surface of the particle, while the velocity component field normal to S
b and the tangential stresses are assumed to vanish on S
b; that is,
(28)
n
^
·
v
=
0
,
r
∈
S
b
,
(29)
n
^
·
Π
~
·
I
~
-
n
^
⊗
n
^
=
0
,
r
∈
S
b
,
where n
^ is the outer unit normal vector. Equations (23)–(29) define a well-posed Happel-type boundary value problem for 3D domains V
R
3, bounded in our case by two arbitrary surfaces S
a and S
b.Papkovich (1932) and Neuber (1934) proposed a differential representation of the flow fields in terms of harmonic functions [7, 8], which is derivable from the well-known Naghdi-Hsu solution [9] and it is applicable to axisymmetric but also to nonaxisymmetric domains as in our case. We profit by the major advantage of the Papkovich-Neuber differential representation, which can be used to obtain solutions of creeping flow in cell models where the shape of the particles is genuinely three-dimensional. Let us notice that the loss of symmetry is caused by the imposed rotation of the particles. Thus, in terms of two harmonic potentials Φ and Ψ,
(30)
Δ
Φ
=
0
,
Δ
Ψ
=
0
,
the Papkovich-Neuber general solution represents the velocity and the total pressure fields appearing in (23) via the action of differential operators on Φ and Ψ, that is,
(31)
v
=
Φ
-
1
2
∇
r
·
Φ
+
Ψ
,
(32)
P
=
P
0
-
μ
f
∇
·
Φ
,
whereas P
0 is a constant pressure of reference usually assigned at a convenient point. It can be easily confirmed that (31) and (32) satisfy Stokes equations (23).Ellipsoidal geometry provides us with the most widely used framework for representing small particles of arbitrary shape embedded within a fluid that flows according to Stokes law. This nonaxisymmetric flow is governed by the genuine 3D ellipsoidal geometry, which embodies the complete anisotropy of the three-dimensional space. The above Happel-type problem (23)–(26), accompanied by the appropriate boundary conditions (27)–(29), is solved with the aim of the Papkovich-Neuber differential representation (30)–(32), using ellipsoidal coordinates. The fields are provided in a closed form fashion as full series expansions of ellipsoidal harmonic eigenfunctions. The velocity, to the first degree, which represents the leading term of the series, is sufficient for most engineering applications and also provides us with the corresponding full 3D solution for the sphere [10] after a proper reduction. The whole analysis is based on the Lamé functions and the theory of ellipsoidal harmonics. In fact, only harmonics of degree less than or equal to two are needed to obtain the velocity field of the first degree. Analytical expressions for the leading terms of the total pressure, the angular velocity, and the stress tensor fields are provided in [5] as well.Introducing the two boundary ellipsoidal surfaces with respect to the radial ellipsoidal variableρ
∈
h
2
,
+
∞ as ρ
=
ρ
a at S
a for the inner on the particle and ρ
=
ρ
b at S
b for the outer boundary of the fluid envelope, then V
R
3
=
ρ
,
μ
,
ν
:
ρ
∈
ρ
a
,
ρ
b
,
μ
∈
h
3
,
h
2
,
ν
∈
0
,
h
3 confines the actual domain of the flow study. In view of the orthonormal unit vector ρ
^ in ellipsoidal coordinates and in terms of the connection formulae between ellipsoidal and spherical harmonics, the mixed Cartesian-ellipsoidal form of the main results obtained in [5] admits more tractable expressions now via the spherical harmonics, where any numerical treatment is much more feasible.After some extended algebra the velocity field renders(33)
v
=
U
+
Ω
×
r
+
Z
ρ
+
∑
κ
=
1
3
‍
H
κ
ρ
E
1
κ
(
r
)
+
ρ
^
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
∑
j
=
1
3
‍
Θ
j
ρ
E
1
j
(
r
)
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
2
i
+
j
(
r
)
,
provided that
(34)
Z
ρ
=
-
E
3
7
ρ
b
U
·
∑
j
=
1
3
‍
I
0
1
ρ
-
I
0
1
ρ
a
+
E
1
j
ρ
a
2
I
1
j
ρ
-
I
1
j
ρ
a
×
x
^
j
⊗
x
^
j
N
j
,
(35)
H
κ
ρ
=
3
2
∑
j
=
1
j
≠
κ
3
‍
h
j
e
j
κ
I
1
κ
ρ
-
I
1
κ
ρ
a
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
-
h
κ
e
κ
j
I
1
j
ρ
-
I
1
j
ρ
a
+
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
×
I
2
j
+
κ
ρ
-
I
2
j
+
κ
ρ
a
x
^
j
h
j
,
κ
=
1,2
,
3
,
whereas
(36)
Θ
κ
ρ
=
-
2
h
κ
E
3
7
ρ
b
h
1
h
2
h
3
U
·
x
^
κ
N
κ
1
-
E
1
κ
ρ
a
E
1
κ
ρ
2
,
κ
=
1,2
,
3
,
(37)
Φ
κ
l
ρ
=
3
h
κ
h
1
h
2
h
3
e
κ
l
1
-
E
1
κ
ρ
a
E
1
κ
ρ
2
1
E
1
l
ρ
2
,
κ
,
l
=
1,2
,
3
,
κ
≠
l
.
On the other hand, the constants N
κ
,
κ
=
1,2
,
3, and e
κ
l, κ
,
l
=
1,2
,
3
,
κ
≠
l, present in (34), (36) and (35), (37), respectively, can be found analytically in [5].Once the above constants are calculated, the velocity field is obtained in terms of the applied fieldsU and Ω via (33). The total pressure field assumes the expression
(38)
P
=
P
0
+
μ
f
E
3
7
ρ
ρ
2
-
μ
2
ρ
2
-
ν
2
×
∑
j
=
1
3
‍
Θ
j
ρ
ρ
2
-
ρ
a
2
E
1
j
+
3
h
1
h
2
h
3
∑
i
,
j
=
1
i
≠
j
3
‍
h
i
e
i
j
E
2
i
+
j
ρ
2
E
2
i
+
j
.
Matching results (25) and (33), we arrive at
(39)
ω
=
Ω
+
h
1
h
2
h
3
2
∑
j
=
1
3
‍
x
^
j
×
H
j
ρ
h
j
+
ρ
2
-
h
3
2
ρ
2
-
h
2
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
×
d
Z
ρ
d
ρ
+
∑
j
=
1
3
‍
d
H
j
ρ
d
ρ
E
1
j
-
h
1
h
2
h
3
4
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
×
∑
j
=
1
3
‍
Θ
j
ρ
h
j
x
^
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
1
i
h
j
x
^
j
+
E
1
j
h
i
x
^
i
for the vorticity field.On the other hand, substitution of (38) and (33) with respect to (24) into (26) reveals that the stress tensor yields
(40)
Π
~
=
ρ
f
g
h
-
P
I
~
-
μ
f
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
^
⊗
d
Z
ρ
d
ρ
+
d
Z
ρ
d
ρ
⊗
ρ
^
+
∑
j
=
1
3
‍
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
1
j
×
ρ
^
⊗
d
H
j
ρ
d
ρ
+
d
H
j
ρ
d
ρ
⊗
ρ
^
+
h
1
h
2
h
3
h
j
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
x
^
j
⊗
H
j
ρ
+
H
j
ρ
⊗
x
^
j
-
μ
f
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
⊗
ρ
^
+
E
1
j
h
i
ρ
^
⊗
x
^
i
+
x
^
i
⊗
ρ
^
×
∑
j
=
1
3
‍
d
Θ
j
ρ
d
ρ
E
1
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
d
Φ
i
j
ρ
d
ρ
E
2
i
+
j
+
h
1
h
2
h
3
2
×
∑
j
=
1
3
‍
Θ
j
ρ
h
j
ρ
^
⊗
x
^
j
+
x
^
j
⊗
ρ
^
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
×
E
1
i
h
j
ρ
^
⊗
x
^
j
+
x
^
j
⊗
ρ
^
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
∑
j
=
1
3
‍
+
E
1
j
h
i
ρ
^
⊗
x
^
i
+
x
^
i
⊗
ρ
^
-
μ
f
2
S
~
∑
j
=
1
3
‍
Θ
j
ρ
E
1
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
2
i
+
j
,
where the dyadic S
~ is given as a function of the metric coefficients in ellipsoidal geometry in the form
(41)
S
~
=
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
h
ρ
-
1
ρ
2
-
μ
2
+
1
ρ
2
-
ν
2
ρ
^
⊗
ρ
^
+
μ
^
⊗
μ
^
ρ
2
-
μ
2
+
ν
^
⊗
ν
^
ρ
2
-
ν
2
+
μ
h
μ
ρ
2
-
μ
2
ρ
^
⊗
μ
^
+
μ
^
⊗
ρ
^
+
ν
h
ν
ρ
2
-
ν
2
ρ
^
⊗
ν
^
+
ν
^
⊗
ρ
^
.
With the aim of the much helpful connection formulas derived in this work, we can interpret the ellipsoidal harmonics E
1
κ and E
2
6
-
κ, κ
=
1,2
,
3, in terms of the spherical harmonics.


### 3.2.  Low-Frequency Electromagnetic Scattering by Impenetrable Ellipsoidal Bodies

Practical applications in several physical areas relative to electromagnetism (e.g., geophysics) are often concerned with the problem of identifying and retrieving metallic anomalies and impenetrable obstacles, which are buried under the surface of the conductive Earth. The detailed way where these problems are dealt with can be found analytically in relevant reference [11] of this particular application.In that sense, we consider a solid body of arbitrary shape with impenetrable surfaceS. The perfectly electrically conducting target is embedded in conductive, homogeneous, isotropic, and nonmagnetic medium of conductivity σ and of permeability μ (usually approximated by the permeability of free space μ
0), where the complex-valued wave number at the operation frequency
(42)
k
=
i
ω
μ
σ
=
ω
μ
σ
2
1
+
i
is provided in terms of imaginary unit i and low circular frequency ω, the dielectric permittivity ɛ being ɛ
≪
σ
/
ω. The external three-dimensional space is W
R
3, considered to be smooth and bounded or not bounded, depending on the physical problem. Without loss of generality, the harmonic time-dependence exp
⁡
-
i
ω
t of all field quantities is implied; thus they are written in terms of the position vector r
=
x
1
,
x
2
,
x
3. The metallic object is illuminated by a known magnetic dipole source m, located at a prescribed position r
0 and arbitrarily orientated, which is defined via
(43)
m
=
∑
j
=
1
3
‍
m
j
x
^
j
.
It is operated at the acceptable low-frequency regime and produces three-dimensional incident waves in an arbitrary direction, the so-called electromagnetic incident fields H
in and E
in, which are radiated by the magnetic dipole (43) and are scattered by the solid target, creating the scattered fields H
s and E
s. Therein,
(44)
H
t
=
H
in
+
H
s
,
E
t
=
E
in
+
E
s
for every r
∈
W
R
3
-
r
0 are the total magnetic and electric fields, given by superposition of corresponding incident and scattered fields, where the singular point r
0 has been excluded.In order to put those tools together within the framework of the low-frequency diffusive scattering theory method [12, 13], we expand the incident (in), the scattered (s), and, consequently, the total (t) electromagnetic fields in a Rayleigh-like manner of positive integral powers of i
k, where k is given by (42), such as
(45)
H
x
=
∑
l
=
0
∞
‍
H
l
x
i
k
l
,
E
x
=
∑
l
=
0
∞
‍
E
l
x
i
k
l
,
x
=
in
,
s
,
t
.
So, Maxwell’s equations
(46)
∇
×
E
x
=
i
ω
μ
H
x
,
∇
×
H
x
=
-
i
ω
ɛ
+
σ
E
x
≅
ɛ
≪
σ
/
ω
σ
E
x
,
x
=
in
,
s
,
t
,
are reduced into low-frequency forms
(47)
σ
∇
×
E
l
x
=
-
H
l
-
2
x
for
l
≥
2
,
∇
×
H
l
x
=
σ
E
l
x
for
l
≥
0
,
x
=
in
,
s
,
t
,
where, either in (46) or in (47), magnetic and electric fields are divergence-free for any r
∈
W
R
3
-
r
0; that is,
(48)
∇
·
H
x
=
∇
·
E
x
=
0
,
∇
·
H
l
x
=
∇
·
E
l
x
=
0
for
l
≥
0
,
x
=
in
,
s
,
t
.
The gradient operator ∇ involved in the above relations operates on r.LettingR
=
r
-
r
0 the electromagnetic incident fields generated by the magnetic dipole m take the expressions [12]
(49)
H
in
=
1
4
π
k
2
+
i
k
R
-
1
R
2
m
-
k
2
+
3
i
k
R
-
3
R
2
r
⊗
r
·
m
R
2
e
i
k
R
R
,
E
in
=
ω
μ
k
4
π
1
+
i
k
R
m
×
r
R
e
i
k
R
R
,
where the symbol ⊗ denotes juxtaposition. Extended algebraic calculations upon (49) lead to the low-frequency forms (45) for the incident fields (x
=
in), where, by implying the hypothesis of low frequency, we conclude that the first four terms of the expansions are adequate (see [11] for details).Therefore, our analysis has confined these important terms of the expansions for the scattered fields, as well. Those are the static term (Rayleigh approximation) forl
=
0 and the dynamic ones for l
=
1,2
,
3, while the terms for l
≥
4 are being expected to be very small and, consequently, they are neglected. Hence, the scattered magnetic field
(50)
H
s
=
H
0
s
+
H
2
s
i
k
2
+
H
3
s
i
k
3
+
O
i
k
4
and the scattered electric field
(51)
E
s
=
E
2
s
i
k
2
+
O
i
k
4
inherit similar forms to those of the incident fields [11], where H
0
s
,
H
2
s
,
H
3
s, and E
2
s are to be evaluated. Substituting the wave number k of the surrounding medium from (42) into relations (50) and (51), trivial analysis yields
(52)
H
s
=
H
0
s
+
ω
μ
σ
ω
μ
σ
2
H
3
s
+
i
ω
μ
σ
ω
μ
σ
2
H
3
s
-
H
2
s
+
O
i
k
4
,
(53)
E
s
=
-
i
ω
μ
σ
E
2
s
+
O
i
k
4
,
respectively. The electric field (53) is purely imaginary-valued, needing only E
2
s, while the magnetic field is complex-valued, noticing that the electromagnetic fields H
2
s and H
0
s are adequate for the full solution, since the contribution of H
3
s, as the outcome of the corresponding constant field, stands for a very small correction to both real and imaginary parts of the scattered magnetic field (52). On the other hand, H
1
s
=
E
0
s
=
E
1
s
=
E
3
s
=
0, in absence of the corresponding incident fields.Straightforward calculations on Maxwell’s equations (47) for x
=
s and elaborate use of identity ∇
×
∇
×
f
=
∇
∇
·
f
-
Δ
f, f being any vector, result in mixed boundary value problems, which are possibly coupled to one another. Those are given by
(54)
Δ
H
0
s
=
0
⟹
H
0
s
=
∇
Φ
0
s
,
since
∇
·
H
0
s
=
0
,
∇
×
H
0
s
=
0
,
(55)
Δ
H
2
s
=
H
0
s
⟹
H
2
s
=
Φ
2
s
+
1
2
r
Φ
0
s
,
σ
E
2
s
=
∇
×
H
2
s
,
since
∇
·
H
2
s
=
∇
·
E
2
s
=
0
,
(56)
Δ
H
3
s
=
0
⟹
H
3
s
=
∇
Φ
3
s
,
since
∇
·
H
3
s
=
0
,
which are written in terms of the harmonic potentials Φ
0
s
,
Φ
2
s, and Φ
3
s that satisfy
(57)
Δ
Φ
0
s
=
Δ
Φ
3
s
=
0
,
Δ
Φ
2
s
=
0
.
It is worth mentioning that for l
=
0,3 standard Laplace equations must be solved for H
0
s and H
3
s fields, whilst the inhomogeneous vector Laplace equation (55), coupled with the solution of (54), is a Poisson partial differential equation; provided that the zero-order scattered field H
0
s is obtained, the second-order scattered field H
2
s can be written as a general vector harmonic function Φ
2
s plus a particular solution (
1
/
2
)
r
Φ
0
s. As for the scattered electric field E
2
s for l
=
2, it is given by the curl of the corresponding magnetic field via (55).The set of problems (54)–(57) has to be solved by using the proper perfectly electrically conducting boundary conditions on the surface S of the perfectly electrically conducting body. They concern the total fields (44) at each order l, where, using outward unit normal vector, the normal component of the total magnetic field and the tangential component of the total electric field are canceling. Hence, combining (44) and (45) we readily obtain
(58)
n
^
·
H
l
in
+
H
l
s
=
0
for
l
=
0,2
,
3
,
n
^
×
E
2
in
+
E
2
s
=
0
,
r
∈
S
,
where the Silver-Müller radiation conditions [12] at infinity must automatically be satisfied.The most arbitrary shape, which is consistent with the well-known orthogonal curvilinear systems we know and reflects the random kind of objects that can be found underground, is the ellipsoidal coordinate system. Therefore, we can readily apply the mathematical tools of the ellipsoidal nonaxisymmetric geometry to our scattering model described earlier and after cumbersome yet rigorous calculations in terms of ellipsoidal harmonic eigenfunctions [1], we may provide the scattered fields as infinite series expansions. This work was accomplished in relevant reference [11], where the fields were presented in terms of the theory of ellipsoidal harmonics and a proper reduction process to the spherical model [14] was demonstrated. At this point, we must mention that, in electromagnetic problems such as ours where the source is far away from the object, the low-frequency approximation is attainable and therein the use of the first degrees of the ellipsoidal harmonic eigenfunctions, in our case let us say n
=
0,1
,
2, is more than enough to reach convergence. Thus, the results collected in [11] use this type of harmonics, which was the key to obtaining closed-type forms for the l-order (l
=
0,1
,
2,3) low-frequency electromagnetic fields (52) and (53), in terms of ellipsoidal harmonics.Introducing the boundary ellipsoidal metal surface with respect to the radial ellipsoidal variableρ
∈
h
2
,
+
∞ as ρ
=
ρ
s
=
a
1 at S, a
1 being the major axis of the ellipsoid under consideration, then the region of electromagnetic scattering is W
R
3
=
ρ
,
μ
,
ν
:
ρ
∈
a
1
,
+
∞
,
μ
∈
h
3
,
h
2
,
ν
∈
0
,
h
3. In view of the orthonormal unit vector ρ
^ in ellipsoidal coordinates and in terms of the connection formulae between ellipsoidal and spherical harmonics, the mixed Cartesian-ellipsoidal form of the main results obtained in [11] provides us with more tractable expressions via the known spherical harmonics, where any numerical interpretation is much easier.Collecting everything together, the scattered magnetic field is(59)
H
s
=
∑
n
=
0
2
‍
∑
m
=
1
2
n
+
1
‍
2
n
+
1
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
b
n
m
+
ω
μ
σ
ω
μ
σ
2
a
n
m
+
i
ω
μ
σ
×
ω
μ
σ
2
a
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
c
n
m
+
b
n
m
2
r
I
n
m
ρ
E
n
m
+
O
i
k
4
,
while for the scattered electric field we have
(60)
E
s
=
-
i
ω
μ
×
∑
n
=
0
2
‍
∑
m
=
1
2
n
+
1
‍
2
n
+
1
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
c
n
m
+
b
n
m
2
r
+
O
i
k
4
.
The coefficients a
n
m, b
n
m, and c
n
m are found explicitly in [11].Moreover, the gradient present in (59) and (60) is easily evaluated since
(61)
∇
r
n
Y
n
m
r
^
=
r
n
n
Y
n
m
(
r
^
)
r
r
^
+
1
sin
2
θ
n
j
n
+
1
m
Y
n
+
1
m
(
r
^
)
-
(
n
+
1
)
j
n
m
Y
n
-
1
m
(
r
^
)
θ
^
Y
n
m
(
r
^
)
r
+
i
m
r
sin
θ
Y
n
m
(
r
^
)
ϕ
^
,
where
(62)
j
n
m
=
n
2
-
m
2
4
n
2
-
1
.
The ellipsoidal harmonics E
n
m comprise the amenable derived spherical harmonics.


## 3.1. Particle-in-Cell Models for Stokes Flow in Ellipsoidal Geometry

A dimensionless criterion, which determines the relative importance of inertial and viscous effects in fluid dynamics, is the Reynolds number (Re). The pseudosteady and nonaxisymmetric creeping flow (Re
≪
1) of an incompressible, viscous fluid with dynamic viscosity μ
f and mass density ρ
f is described by the well-known Stokes equations.Let us consider any smooth and bounded or unbounded three-dimensional domainV
R
3, depending on the physical problem, in terms of the position vector r
=
x
1
,
x
2
,
x
3. Therein, Stokes flow governing equations connect the biharmonic vector velocity v with the harmonic scalar total pressure field P via
(23)
μ
f
Δ
v
=
∇
P
,
∇
·
v
=
0
for every r
∈
V
R
3.The left-hand side of (23) states that, in creeping flow, the viscous force compensates the force caused by the pressure gradient, while the right-hand side of it secures the incompressibility of the fluid. The fluid total pressure is related to the thermodynamic pressure p through
(24)
P
=
p
+
ρ
f
g
h
,
where the contribution of the term ρ
f
g
h, with g being the acceleration of gravity, refers to the gravitational pressure force, corresponding to a height of reference h. Once the velocity field is obtained, the harmonic vorticity field ω is defined as
(25)
ω
=
1
2
∇
×
v
,
while the stress dyadic Π
~ is assumed to be
(26)
Π
~
=
-
p
I
~
+
μ
f
∇
⊗
v
+
∇
⊗
v
⊺
in terms of the unit dyadic I
~
=
∑
j
=
1
3
‍
x
^
j
⊗
x
^
j, where the symbols ⊗ and ⊺ denote juxtaposition and transposition, respectively.One of the important areas of applications concerns the construction of particle-in-cell models, which are useful in the development of simple but reliable analytical expressions for heat and mass transfer in swarms of particles in the case of concentrated suspensions. In applied type analysis it is not usually necessary to have detailed solution of the flow field over the entire swarm of particles taking into account the exact positions of the particles, since such solutions are cumbersome to use. Thus, the technique of cell models is adopted where the mathematical treatment of each problem is based on the assumption that a three-dimensional assemblage may be considered to consist of a number of identical unit cells. Each of these cells contains a particle surrounded by a fluid envelope, containing a volume of fluid sufficient to make the fractional void volume in the cell identical to that in the entire assemblage.With the aim of mathematical modeling of the particular physical problem, we inherit, chosen among others cited in original paper [5] from where the main results were drawn, the Happel cell model [6] in which both the particle and the outer envelope enjoy spherical symmetry. In view of the Happel-type model [5, 6], we consider a fluid-particle system consisting of any finite number of rigid particles of arbitrary shape. Introducing the particle-in-cell model, we examine the Stokes flow of one of the assemblages of particles neglecting the interaction with other particles or with the bounded walls of a container. Let S
a denote the surface of the particle of the swarm, which is solid, is moving with a known constant translational velocity U
=
(
U
1
,
U
2
,
U
3
) in an arbitrary direction, and is rotating, also arbitrarily, with a defined constant angular velocity Ω
=
(
Ω
1
,
Ω
2
,
Ω
3
). It lives within an otherwise quiescent fluid layer, which is confined by the outer surface denoted by S
b. Following the formulation of Happel [6], extended to three-dimensional flows, the necessary nonslip flow conditions
(27)
v
=
U
+
Ω
×
r
,
r
∈
S
a
are imposed on the surface of the particle, while the velocity component field normal to S
b and the tangential stresses are assumed to vanish on S
b; that is,
(28)
n
^
·
v
=
0
,
r
∈
S
b
,
(29)
n
^
·
Π
~
·
I
~
-
n
^
⊗
n
^
=
0
,
r
∈
S
b
,
where n
^ is the outer unit normal vector. Equations (23)–(29) define a well-posed Happel-type boundary value problem for 3D domains V
R
3, bounded in our case by two arbitrary surfaces S
a and S
b.Papkovich (1932) and Neuber (1934) proposed a differential representation of the flow fields in terms of harmonic functions [7, 8], which is derivable from the well-known Naghdi-Hsu solution [9] and it is applicable to axisymmetric but also to nonaxisymmetric domains as in our case. We profit by the major advantage of the Papkovich-Neuber differential representation, which can be used to obtain solutions of creeping flow in cell models where the shape of the particles is genuinely three-dimensional. Let us notice that the loss of symmetry is caused by the imposed rotation of the particles. Thus, in terms of two harmonic potentials Φ and Ψ,
(30)
Δ
Φ
=
0
,
Δ
Ψ
=
0
,
the Papkovich-Neuber general solution represents the velocity and the total pressure fields appearing in (23) via the action of differential operators on Φ and Ψ, that is,
(31)
v
=
Φ
-
1
2
∇
r
·
Φ
+
Ψ
,
(32)
P
=
P
0
-
μ
f
∇
·
Φ
,
whereas P
0 is a constant pressure of reference usually assigned at a convenient point. It can be easily confirmed that (31) and (32) satisfy Stokes equations (23).Ellipsoidal geometry provides us with the most widely used framework for representing small particles of arbitrary shape embedded within a fluid that flows according to Stokes law. This nonaxisymmetric flow is governed by the genuine 3D ellipsoidal geometry, which embodies the complete anisotropy of the three-dimensional space. The above Happel-type problem (23)–(26), accompanied by the appropriate boundary conditions (27)–(29), is solved with the aim of the Papkovich-Neuber differential representation (30)–(32), using ellipsoidal coordinates. The fields are provided in a closed form fashion as full series expansions of ellipsoidal harmonic eigenfunctions. The velocity, to the first degree, which represents the leading term of the series, is sufficient for most engineering applications and also provides us with the corresponding full 3D solution for the sphere [10] after a proper reduction. The whole analysis is based on the Lamé functions and the theory of ellipsoidal harmonics. In fact, only harmonics of degree less than or equal to two are needed to obtain the velocity field of the first degree. Analytical expressions for the leading terms of the total pressure, the angular velocity, and the stress tensor fields are provided in [5] as well.Introducing the two boundary ellipsoidal surfaces with respect to the radial ellipsoidal variableρ
∈
h
2
,
+
∞ as ρ
=
ρ
a at S
a for the inner on the particle and ρ
=
ρ
b at S
b for the outer boundary of the fluid envelope, then V
R
3
=
ρ
,
μ
,
ν
:
ρ
∈
ρ
a
,
ρ
b
,
μ
∈
h
3
,
h
2
,
ν
∈
0
,
h
3 confines the actual domain of the flow study. In view of the orthonormal unit vector ρ
^ in ellipsoidal coordinates and in terms of the connection formulae between ellipsoidal and spherical harmonics, the mixed Cartesian-ellipsoidal form of the main results obtained in [5] admits more tractable expressions now via the spherical harmonics, where any numerical treatment is much more feasible.After some extended algebra the velocity field renders(33)
v
=
U
+
Ω
×
r
+
Z
ρ
+
∑
κ
=
1
3
‍
H
κ
ρ
E
1
κ
(
r
)
+
ρ
^
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
∑
j
=
1
3
‍
Θ
j
ρ
E
1
j
(
r
)
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
2
i
+
j
(
r
)
,
provided that
(34)
Z
ρ
=
-
E
3
7
ρ
b
U
·
∑
j
=
1
3
‍
I
0
1
ρ
-
I
0
1
ρ
a
+
E
1
j
ρ
a
2
I
1
j
ρ
-
I
1
j
ρ
a
×
x
^
j
⊗
x
^
j
N
j
,
(35)
H
κ
ρ
=
3
2
∑
j
=
1
j
≠
κ
3
‍
h
j
e
j
κ
I
1
κ
ρ
-
I
1
κ
ρ
a
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
-
h
κ
e
κ
j
I
1
j
ρ
-
I
1
j
ρ
a
+
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
h
j
e
j
κ
E
1
j
ρ
a
2
+
h
κ
e
κ
j
E
1
κ
ρ
a
2
×
I
2
j
+
κ
ρ
-
I
2
j
+
κ
ρ
a
x
^
j
h
j
,
κ
=
1,2
,
3
,
whereas
(36)
Θ
κ
ρ
=
-
2
h
κ
E
3
7
ρ
b
h
1
h
2
h
3
U
·
x
^
κ
N
κ
1
-
E
1
κ
ρ
a
E
1
κ
ρ
2
,
κ
=
1,2
,
3
,
(37)
Φ
κ
l
ρ
=
3
h
κ
h
1
h
2
h
3
e
κ
l
1
-
E
1
κ
ρ
a
E
1
κ
ρ
2
1
E
1
l
ρ
2
,
κ
,
l
=
1,2
,
3
,
κ
≠
l
.
On the other hand, the constants N
κ
,
κ
=
1,2
,
3, and e
κ
l, κ
,
l
=
1,2
,
3
,
κ
≠
l, present in (34), (36) and (35), (37), respectively, can be found analytically in [5].Once the above constants are calculated, the velocity field is obtained in terms of the applied fieldsU and Ω via (33). The total pressure field assumes the expression
(38)
P
=
P
0
+
μ
f
E
3
7
ρ
ρ
2
-
μ
2
ρ
2
-
ν
2
×
∑
j
=
1
3
‍
Θ
j
ρ
ρ
2
-
ρ
a
2
E
1
j
+
3
h
1
h
2
h
3
∑
i
,
j
=
1
i
≠
j
3
‍
h
i
e
i
j
E
2
i
+
j
ρ
2
E
2
i
+
j
.
Matching results (25) and (33), we arrive at
(39)
ω
=
Ω
+
h
1
h
2
h
3
2
∑
j
=
1
3
‍
x
^
j
×
H
j
ρ
h
j
+
ρ
2
-
h
3
2
ρ
2
-
h
2
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
×
d
Z
ρ
d
ρ
+
∑
j
=
1
3
‍
d
H
j
ρ
d
ρ
E
1
j
-
h
1
h
2
h
3
4
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
×
∑
j
=
1
3
‍
Θ
j
ρ
h
j
x
^
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
1
i
h
j
x
^
j
+
E
1
j
h
i
x
^
i
for the vorticity field.On the other hand, substitution of (38) and (33) with respect to (24) into (26) reveals that the stress tensor yields
(40)
Π
~
=
ρ
f
g
h
-
P
I
~
-
μ
f
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
^
⊗
d
Z
ρ
d
ρ
+
d
Z
ρ
d
ρ
⊗
ρ
^
+
∑
j
=
1
3
‍
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
1
j
×
ρ
^
⊗
d
H
j
ρ
d
ρ
+
d
H
j
ρ
d
ρ
⊗
ρ
^
+
h
1
h
2
h
3
h
j
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
x
^
j
⊗
H
j
ρ
+
H
j
ρ
⊗
x
^
j
-
μ
f
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
ρ
^
⊗
ρ
^
+
E
1
j
h
i
ρ
^
⊗
x
^
i
+
x
^
i
⊗
ρ
^
×
∑
j
=
1
3
‍
d
Θ
j
ρ
d
ρ
E
1
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
d
Φ
i
j
ρ
d
ρ
E
2
i
+
j
+
h
1
h
2
h
3
2
×
∑
j
=
1
3
‍
Θ
j
ρ
h
j
ρ
^
⊗
x
^
j
+
x
^
j
⊗
ρ
^
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
×
E
1
i
h
j
ρ
^
⊗
x
^
j
+
x
^
j
⊗
ρ
^
ρ
2
-
h
3
2
ρ
2
-
h
2
2
ρ
2
-
μ
2
ρ
2
-
ν
2
∑
j
=
1
3
‍
+
E
1
j
h
i
ρ
^
⊗
x
^
i
+
x
^
i
⊗
ρ
^
-
μ
f
2
S
~
∑
j
=
1
3
‍
Θ
j
ρ
E
1
j
+
∑
i
,
j
=
1
i
≠
j
3
‍
Φ
i
j
ρ
E
2
i
+
j
,
where the dyadic S
~ is given as a function of the metric coefficients in ellipsoidal geometry in the form
(41)
S
~
=
2
ρ
2
-
μ
2
ρ
2
-
ν
2
×
ρ
h
ρ
-
1
ρ
2
-
μ
2
+
1
ρ
2
-
ν
2
ρ
^
⊗
ρ
^
+
μ
^
⊗
μ
^
ρ
2
-
μ
2
+
ν
^
⊗
ν
^
ρ
2
-
ν
2
+
μ
h
μ
ρ
2
-
μ
2
ρ
^
⊗
μ
^
+
μ
^
⊗
ρ
^
+
ν
h
ν
ρ
2
-
ν
2
ρ
^
⊗
ν
^
+
ν
^
⊗
ρ
^
.
With the aim of the much helpful connection formulas derived in this work, we can interpret the ellipsoidal harmonics E
1
κ and E
2
6
-
κ, κ
=
1,2
,
3, in terms of the spherical harmonics.


## 3.2.  Low-Frequency Electromagnetic Scattering by Impenetrable Ellipsoidal Bodies

Practical applications in several physical areas relative to electromagnetism (e.g., geophysics) are often concerned with the problem of identifying and retrieving metallic anomalies and impenetrable obstacles, which are buried under the surface of the conductive Earth. The detailed way where these problems are dealt with can be found analytically in relevant reference [11] of this particular application.In that sense, we consider a solid body of arbitrary shape with impenetrable surfaceS. The perfectly electrically conducting target is embedded in conductive, homogeneous, isotropic, and nonmagnetic medium of conductivity σ and of permeability μ (usually approximated by the permeability of free space μ
0), where the complex-valued wave number at the operation frequency
(42)
k
=
i
ω
μ
σ
=
ω
μ
σ
2
1
+
i
is provided in terms of imaginary unit i and low circular frequency ω, the dielectric permittivity ɛ being ɛ
≪
σ
/
ω. The external three-dimensional space is W
R
3, considered to be smooth and bounded or not bounded, depending on the physical problem. Without loss of generality, the harmonic time-dependence exp
⁡
-
i
ω
t of all field quantities is implied; thus they are written in terms of the position vector r
=
x
1
,
x
2
,
x
3. The metallic object is illuminated by a known magnetic dipole source m, located at a prescribed position r
0 and arbitrarily orientated, which is defined via
(43)
m
=
∑
j
=
1
3
‍
m
j
x
^
j
.
It is operated at the acceptable low-frequency regime and produces three-dimensional incident waves in an arbitrary direction, the so-called electromagnetic incident fields H
in and E
in, which are radiated by the magnetic dipole (43) and are scattered by the solid target, creating the scattered fields H
s and E
s. Therein,
(44)
H
t
=
H
in
+
H
s
,
E
t
=
E
in
+
E
s
for every r
∈
W
R
3
-
r
0 are the total magnetic and electric fields, given by superposition of corresponding incident and scattered fields, where the singular point r
0 has been excluded.In order to put those tools together within the framework of the low-frequency diffusive scattering theory method [12, 13], we expand the incident (in), the scattered (s), and, consequently, the total (t) electromagnetic fields in a Rayleigh-like manner of positive integral powers of i
k, where k is given by (42), such as
(45)
H
x
=
∑
l
=
0
∞
‍
H
l
x
i
k
l
,
E
x
=
∑
l
=
0
∞
‍
E
l
x
i
k
l
,
x
=
in
,
s
,
t
.
So, Maxwell’s equations
(46)
∇
×
E
x
=
i
ω
μ
H
x
,
∇
×
H
x
=
-
i
ω
ɛ
+
σ
E
x
≅
ɛ
≪
σ
/
ω
σ
E
x
,
x
=
in
,
s
,
t
,
are reduced into low-frequency forms
(47)
σ
∇
×
E
l
x
=
-
H
l
-
2
x
for
l
≥
2
,
∇
×
H
l
x
=
σ
E
l
x
for
l
≥
0
,
x
=
in
,
s
,
t
,
where, either in (46) or in (47), magnetic and electric fields are divergence-free for any r
∈
W
R
3
-
r
0; that is,
(48)
∇
·
H
x
=
∇
·
E
x
=
0
,
∇
·
H
l
x
=
∇
·
E
l
x
=
0
for
l
≥
0
,
x
=
in
,
s
,
t
.
The gradient operator ∇ involved in the above relations operates on r.LettingR
=
r
-
r
0 the electromagnetic incident fields generated by the magnetic dipole m take the expressions [12]
(49)
H
in
=
1
4
π
k
2
+
i
k
R
-
1
R
2
m
-
k
2
+
3
i
k
R
-
3
R
2
r
⊗
r
·
m
R
2
e
i
k
R
R
,
E
in
=
ω
μ
k
4
π
1
+
i
k
R
m
×
r
R
e
i
k
R
R
,
where the symbol ⊗ denotes juxtaposition. Extended algebraic calculations upon (49) lead to the low-frequency forms (45) for the incident fields (x
=
in), where, by implying the hypothesis of low frequency, we conclude that the first four terms of the expansions are adequate (see [11] for details).Therefore, our analysis has confined these important terms of the expansions for the scattered fields, as well. Those are the static term (Rayleigh approximation) forl
=
0 and the dynamic ones for l
=
1,2
,
3, while the terms for l
≥
4 are being expected to be very small and, consequently, they are neglected. Hence, the scattered magnetic field
(50)
H
s
=
H
0
s
+
H
2
s
i
k
2
+
H
3
s
i
k
3
+
O
i
k
4
and the scattered electric field
(51)
E
s
=
E
2
s
i
k
2
+
O
i
k
4
inherit similar forms to those of the incident fields [11], where H
0
s
,
H
2
s
,
H
3
s, and E
2
s are to be evaluated. Substituting the wave number k of the surrounding medium from (42) into relations (50) and (51), trivial analysis yields
(52)
H
s
=
H
0
s
+
ω
μ
σ
ω
μ
σ
2
H
3
s
+
i
ω
μ
σ
ω
μ
σ
2
H
3
s
-
H
2
s
+
O
i
k
4
,
(53)
E
s
=
-
i
ω
μ
σ
E
2
s
+
O
i
k
4
,
respectively. The electric field (53) is purely imaginary-valued, needing only E
2
s, while the magnetic field is complex-valued, noticing that the electromagnetic fields H
2
s and H
0
s are adequate for the full solution, since the contribution of H
3
s, as the outcome of the corresponding constant field, stands for a very small correction to both real and imaginary parts of the scattered magnetic field (52). On the other hand, H
1
s
=
E
0
s
=
E
1
s
=
E
3
s
=
0, in absence of the corresponding incident fields.Straightforward calculations on Maxwell’s equations (47) for x
=
s and elaborate use of identity ∇
×
∇
×
f
=
∇
∇
·
f
-
Δ
f, f being any vector, result in mixed boundary value problems, which are possibly coupled to one another. Those are given by
(54)
Δ
H
0
s
=
0
⟹
H
0
s
=
∇
Φ
0
s
,
since
∇
·
H
0
s
=
0
,
∇
×
H
0
s
=
0
,
(55)
Δ
H
2
s
=
H
0
s
⟹
H
2
s
=
Φ
2
s
+
1
2
r
Φ
0
s
,
σ
E
2
s
=
∇
×
H
2
s
,
since
∇
·
H
2
s
=
∇
·
E
2
s
=
0
,
(56)
Δ
H
3
s
=
0
⟹
H
3
s
=
∇
Φ
3
s
,
since
∇
·
H
3
s
=
0
,
which are written in terms of the harmonic potentials Φ
0
s
,
Φ
2
s, and Φ
3
s that satisfy
(57)
Δ
Φ
0
s
=
Δ
Φ
3
s
=
0
,
Δ
Φ
2
s
=
0
.
It is worth mentioning that for l
=
0,3 standard Laplace equations must be solved for H
0
s and H
3
s fields, whilst the inhomogeneous vector Laplace equation (55), coupled with the solution of (54), is a Poisson partial differential equation; provided that the zero-order scattered field H
0
s is obtained, the second-order scattered field H
2
s can be written as a general vector harmonic function Φ
2
s plus a particular solution (
1
/
2
)
r
Φ
0
s. As for the scattered electric field E
2
s for l
=
2, it is given by the curl of the corresponding magnetic field via (55).The set of problems (54)–(57) has to be solved by using the proper perfectly electrically conducting boundary conditions on the surface S of the perfectly electrically conducting body. They concern the total fields (44) at each order l, where, using outward unit normal vector, the normal component of the total magnetic field and the tangential component of the total electric field are canceling. Hence, combining (44) and (45) we readily obtain
(58)
n
^
·
H
l
in
+
H
l
s
=
0
for
l
=
0,2
,
3
,
n
^
×
E
2
in
+
E
2
s
=
0
,
r
∈
S
,
where the Silver-Müller radiation conditions [12] at infinity must automatically be satisfied.The most arbitrary shape, which is consistent with the well-known orthogonal curvilinear systems we know and reflects the random kind of objects that can be found underground, is the ellipsoidal coordinate system. Therefore, we can readily apply the mathematical tools of the ellipsoidal nonaxisymmetric geometry to our scattering model described earlier and after cumbersome yet rigorous calculations in terms of ellipsoidal harmonic eigenfunctions [1], we may provide the scattered fields as infinite series expansions. This work was accomplished in relevant reference [11], where the fields were presented in terms of the theory of ellipsoidal harmonics and a proper reduction process to the spherical model [14] was demonstrated. At this point, we must mention that, in electromagnetic problems such as ours where the source is far away from the object, the low-frequency approximation is attainable and therein the use of the first degrees of the ellipsoidal harmonic eigenfunctions, in our case let us say n
=
0,1
,
2, is more than enough to reach convergence. Thus, the results collected in [11] use this type of harmonics, which was the key to obtaining closed-type forms for the l-order (l
=
0,1
,
2,3) low-frequency electromagnetic fields (52) and (53), in terms of ellipsoidal harmonics.Introducing the boundary ellipsoidal metal surface with respect to the radial ellipsoidal variableρ
∈
h
2
,
+
∞ as ρ
=
ρ
s
=
a
1 at S, a
1 being the major axis of the ellipsoid under consideration, then the region of electromagnetic scattering is W
R
3
=
ρ
,
μ
,
ν
:
ρ
∈
a
1
,
+
∞
,
μ
∈
h
3
,
h
2
,
ν
∈
0
,
h
3. In view of the orthonormal unit vector ρ
^ in ellipsoidal coordinates and in terms of the connection formulae between ellipsoidal and spherical harmonics, the mixed Cartesian-ellipsoidal form of the main results obtained in [11] provides us with more tractable expressions via the known spherical harmonics, where any numerical interpretation is much easier.Collecting everything together, the scattered magnetic field is(59)
H
s
=
∑
n
=
0
2
‍
∑
m
=
1
2
n
+
1
‍
2
n
+
1
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
b
n
m
+
ω
μ
σ
ω
μ
σ
2
a
n
m
+
i
ω
μ
σ
×
ω
μ
σ
2
a
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
c
n
m
+
b
n
m
2
r
I
n
m
ρ
E
n
m
+
O
i
k
4
,
while for the scattered electric field we have
(60)
E
s
=
-
i
ω
μ
×
∑
n
=
0
2
‍
∑
m
=
1
2
n
+
1
‍
2
n
+
1
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
I
n
m
ρ
∇
E
n
m
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
-
ρ
^
E
n
m
ρ
2
ρ
2
-
μ
2
ρ
2
-
ν
2
E
n
m
×
c
n
m
+
b
n
m
2
r
+
O
i
k
4
.
The coefficients a
n
m, b
n
m, and c
n
m are found explicitly in [11].Moreover, the gradient present in (59) and (60) is easily evaluated since
(61)
∇
r
n
Y
n
m
r
^
=
r
n
n
Y
n
m
(
r
^
)
r
r
^
+
1
sin
2
θ
n
j
n
+
1
m
Y
n
+
1
m
(
r
^
)
-
(
n
+
1
)
j
n
m
Y
n
-
1
m
(
r
^
)
θ
^
Y
n
m
(
r
^
)
r
+
i
m
r
sin
θ
Y
n
m
(
r
^
)
ϕ
^
,
where
(62)
j
n
m
=
n
2
-
m
2
4
n
2
-
1
.
The ellipsoidal harmonics E
n
m comprise the amenable derived spherical harmonics.


## 4. Conclusions

We introduce explicit formulas relating ellipsoidal harmonics in terms of spherical harmonics up to the third degree, namely, sixteen harmonics in total, valid in the interior and exterior of a confocal triaxial ellipsoid. For a number of physical problems, depending on the prescribed conditions, the existence of source terms, et cetera, these sixteen ellipsoidal harmonics are adequate to ensure convergence of the related expansion to an extent reaching almost over 95 per cent. For these kinds of problems (two applications are advertised in Section3) the expressions supplied can be used directly without prior knowledge of the theory of ellipsoidal harmonics.For the sake of completeness, we indicate that an essentially identical process can be commenced correlating spherical harmonics with the associated ellipsoidal harmonics, which, however, in view of practical problems is meaningless.



---

*Source: 572458-2015-02-02.xml*