Multipole matrix elements of Green function of Laplace equation

0.13 MB·English

by  Karol Makuch

Multipole matrix elements of Green function of Laplace equation
Karol Makuch
Institute of Theoretical Physics, Faculty of Physics,
University of Warsaw, Pasteura 5, 02-093 Warsaw, Poland∗
Przemysław Górka
Department of Mathematics and Information Sciences,
Warsaw University of Technology, Koszykowa 75, 00-662 Warsaw, Poland
(Date textdate;Received textdate;Revised textdate;Accepted textdate;Published textdate)
Multipole matrix elements of Green function of Laplace equation are calculated. The multipole
matrixelementsofGreenfunctioninelectrostaticsdescribepotentialonaspherewhichisproduced
5 byachargedistributedonthesurfaceofadifferent(possiblyoverlapping)sphereofthesameradius.
1 The matrix elements are defined by double convolution of two spherical harmonics with the Green
0 function of Laplace equation. The method we use relies on the fact that in the Fourier space the
2 double convolution has simple form. Therefore we calculate the multipole matrix from its Fourier
n transform. AnimportantpartofourconsiderationsissimplificationofthethreedimensionalFourier
a transformationofgeneralmultipolematrixbyitsrotationalsymmetrytotheone-dimensionalHankel
J transformation.

2

]
h
p
-
h
t
a
m
[

1
v
6
9
3
0
0
.
1
0
5
1
:
v
i
X
r
a
∗ [email protected]
2
I. INTRODUCTION
Metallic spheres in a vacuum, minute spherical particles in a fluid, and spherical inclusions in a solid body are
often considered in statistical physics of dispersive media [1]. Equations which describe those systems are linear and
have spherical symmetry. It is the case of Laplace equation in dielectrics [2], Stokes equations in suspensions [3] and
Laméequations [4]describingsolidbody. The firststepinconsiderationsofdispersivemedia isoften tosolvea single
particle problem. In dielectrics it means to find distribution of charge on the surface of a single metallic sphere in
an external electrostatic potential. To find the solution of single particle problem it is convenient to take spherical
symmetryinto account[5]. Tothis endone introducesasetofscalarfunctions onthe spherewhichis invariantunder
rotations. That is how spherical harmonics Y (ˆr) enter calculations.
lm
Topassfromconsiderationsofasingleparticletothe analysisofmanyparticlesitdemandstoanswerthefollowing
question. What is the electrostatic potential produced by a charge distributed in one spherical surface in the area
occupiedby a differentsphere. The answerto the abovequestioncanbe inferredfromthe followingmultipole matrix
element
[G(R,R′)]
lm,l′m′
1 1
=ZR3d3rZR3d3r′aδ(|r−R|−a)φ∗lm(r−R)G(r−r′)aδ(|r′−R′|−a)φl′m′(r′−R′). (1)
In the above definition G(r) is the Green function of Laplace equation [2],
1
G(r)= ,
4πr
aisradiusofparticles,andone-dimensionalDiracdeltadistributionδ(x)isused. Moreoverφ (r)aresolidharmonics
lm
defined by
φ (r)=rlY (ˆr),
lm lm
with spherical harmonics Y (ˆr) numbered by order l = 0,1,... and azimuthal number m = l,...,l [6]. In our
lm
notation an argument of spherical harmonics is a versor ˆr(θ,φ) in direction described by ang−les θ,φ in spherical
coordinates. Dirac delta distributions in equation (1) reduce three-dimensional integrations to integrals over the
surface of the spheres with radius a centered at positions R and R′.
We can differentiate between two situations. The first situation corresponds to nonoverlapping spheres, i.e. when
R R′ > 2a. In this case the matrix elements defined by equation (1) can be inferred from the results in the
| − |
literature [7]. They have application e.g. in numerical simulations. The second situation is the case of overlapping
configurations, R R′ <2a. Eveniftheparticlesinthesystemcannotoverlaptheremayappearaneedofcalculation
ofoverlappingc|onfi−gura|tionsofthe multipole matrixelements [G(R,R′)] . Forexample,ithasbeenrecognized
lm,l′m′
thatoverlappingconfigurationsappearinmicroscopicexplanationofthefamousClauius-Mossottiformulafordielectric
constant[8]. Inthiscontextintegral d3R′[G(R,R′)]exp ikR′ forlowmultipolenumbersl,l′ has
|R−R′|<2a − lm,l′m′
been considered. Therefore the mulRtipole matrix elements [G(R,R(cid:0)′)] (cid:1) for overlapping configurations play an
lm,l′m′
important role in statistical physics considerations of dispersive media.
In this article we give general expression for the multipole matrix elements [G(R,R′)] defined by equation
lm,l′m′
(1). The new contribution of the current article is the calculation of [G(R,R′)] for overlapping configurations
lm,l′m′
R R′ <2a. We are going to use the result for overlapping configurations in our statistical physics considerations
| − |
of transport properties of dispersive media in further work.
The method of calculation of the multipole matrix elements defined by equation (1) is the following. We observe
that [G(R,R′)] has the form of a double convolution of three functions - two solid harmonics and the Green
lm,l′m′
function. ThereforewecalculatetheFouriertransformofthethreefunctions,taketheirproduct,andthenperformthe
inverseFouriertransformto obtainG(R,R′). An importantelementofourcalculationsisto use sphericalsymmetry
whichallowstoreducethethree-dimensionalFouriertransformofthemultipolematrixtotheone-dimensionalHankel
transform.
II. MULTIPOLE MATRIX IN THE FOURIER SPACE
The aim of this article is calculation of integral (1). We can simplify it using homogeneity of the Laplace equation
which implies that the multipole matrix [G(R,R′)] depends on the relative positions R R′. Therefore from
lm,l′m′ −
now on we consider multipole matrix [G(R)] defined by
lm,l′m′
3
[G(R R′)] =[G(R,R′)]
− lm,l′m′ lm,l′m′
depending on the relative positions.
Next, we observe that the multipole matrix given by the formula (1) has a form of a double convolution of three
ffuunnccttiioonnsG. (Trh).eTthhreeFeofuurniecrtiotrnasnsafroermωl′omf′a(rd)ou≡bleδ(c|orn|v−olau)tφiol′nmo′f(rt)h/eat,hcreoenjfuugnactteiodnfsuinscgtiiovnenωbl∗ymt(hre),parnoddutchteofGtrheeeinr
Fourier transforms, therefore
hGˆ(k)ilm,l′m′ =ωˆl∗m(k)Gˆ (k)ωˆl′m′(k).
In our calculations we use the following definition of the three-dimensional Fourier transform
Gˆlm,l′m′(k)= d3R exp( ikR)Glm,l′m′(R), (2)
Z −
with the inverse transformation given by the formula
1
Glm,l′m′(R)= (2π)3 Z d3k exp(ikR)Gˆlm,l′m′(k).
In this situation the Fourier transforms of ω (r) and G(r) are given respectively by
lm
Gˆ (k)= 1 ,
k2
and
ωˆ (k)=4πal+1( i)lj (ka)Y kˆ .
lm l lm
− (cid:16) (cid:17)
Thelastexpressioncanbecalculatedwiththeuseofequation(5.8.3)fromRef. [6]andwiththeuseoforthonormality
of spherical harmonics. The expression contains spherical Bessel functions j (x) of the order l. Finally, the Fourier
l
transform of multipole matrix (1) is given by the following formula
hGˆ(k)ilm,l′m′ =(4π)2(−i)−l+l′al+l′+2jl(ka)kj2l′(ka)Yl∗m(cid:16)kˆ(cid:17)Yl′m′(cid:16)kˆ(cid:17). (3)
According to the above expression, the Fourier transform of [G(R)] is given by the spherical harmonics and
lm,l′m′
the spherical Bessel functions. At this point we face the main difficulty in our calculations of [G(R)] . The
lm,l′m′
inverse Fourier transform of the above formula has to be calculated. To this end we consider rotational symmetry of
the multipole matrix.
III. ROTATIONAL SYMMETRY OF A MULTIPOLE MATRIX
In the definition of the multipole matrix (1) there appear solid harmonics φ (r) and Green function G(r).
lm
Transformation of solid harmonics under rotation is given by the formula
l
∗
φ (D(α,γ,β)r)= D(l)(α,γ,β) φ (r). (4)
lm lm1
mX1=−lh imm1
The above expression can be inferred from Ref. [6] from which we adopt notation in this article. There, formula
(4.1.1) defines three-dimensional rotation matrix D(α,β,γ) characterized by the Euler angles α,β,γ. Moreover,
D(l) (α,β,γ) denotes the Wigner matrix. Isotropy of the Laplace equation implies that its Green function G(r) is
mm′
invariant under rotation, i.e. G(D(α,β,γ)r)=G(r).
4
Changing the variables in the integrals in Eq. (1) and the above properties of the solid harmonics and the Green
function lead us to the following transformationof the multipole matrix elements under rotation
Glm,l′m′(D(α,β,γ)R)= Dm(l),m1(α,β,γ)(cid:20)Dm(l′′),m′1(α,β,γ)(cid:21)∗Glm1,l′m′1(R). (5)
mX1,m′1
ItisworthmentioningherethattheabovetransformationappliedforRinz direction,R=Rˆe ,andforanyrotations
z
aroundaxis z, thus for β =0 and any other Euler angles, implies that the multipole matrix is diagonalin indexes m,
i.e.
Glm,l′m′(Rˆez)=δm,m′Glm,l′m(Rˆez),
where Kronecker delta δm,m′ appears. It is easy to prove similar diagonality in case of the Fourier transform,
Gˆlm,l′m′(kˆez)=δm,m′Gˆlm,l′m(kˆez), (6)
which appears as a result of Fourier transformationof equation (5) and considerations for proper rotations and wave
vector k=kˆe .
z
IV. FOURIER TRANSFORM OF MULTIPOLE MATRIX - SIMPLIFICATION BY ROTATIONAL
SYMMETRY
ThekeypointofourcalculationsissimplificationoftheFouriertransformofamatrixsatisfyingsymmetryproperty
givenbyequation(5). WesimplifytheFouriertransformforawavevectoralongzdirectionbecausethenthemultipole
matrix is diagonalin indexesm, asis shownby equation(6). For the case k=kˆe expression(2) writteninspherical
z
coordinates reads
∞ π 2π
Gˆlm,l′m′(kˆez)= dR dθ dφR2sinθexp( ikRcosθ)Glm,l′m′(R(R,θ,φ)). (7)
Z Z Z −
0 0 0
We express vector R(R,θ,φ) by a product of the vector Rˆe and the rotation matrix D(α,β,γ) characterized with
z
proper Euler angles. The angles can be deduced from formula (4.1.1) in Ref. [6]. These angles are β = θ, γ = φ,
and any α, e.g. α = 0. Therefore R(R,θ,φ) = D(0, θ, φ)Rˆe . For this rotation and R=Rˆe we us−e symme−try
z z
− −
property (5) which leads to the following expression for multipole matrix G
Glm,l′m′(R(R,θ,φ))= Dm(l),m1(0,−θ,−φ)(cid:20)Dm(l′′),m′1(0,−θ,−φ)(cid:21)∗Glm1,l′m′1(Rˆez). (8)
mX1,m′1
Inthenextstepwerepresentexp( ikRcosθ)fromexpression(7)intheformofaninfiniteseriesofsphericalBessel
−
functions
∞
exp( ikRcosθ)= ( i)l1(2l +1)j (kR)D(l1)(0, θ, φ),number
− − 1 l1 00 − −
lX1=0
whichisdeducedfromformula(5.8.1)and(4.1.26)inthereference[6]. Takingintoconsiderationthelasttwoformulae
in expression (7) yields
∞ ∞
Gˆlm,l′m′(kˆez)=Z dRR2 (−i)l1(2l1+1)jl1(kR)Glm1,l′m′1(Rˆez)×
0 lX1=0
πdθ 2πdφsinθD(l1)(0,θ,φ)D(l) (0,θ,φ) D(l′) (0,θ,φ) ∗. (9)
Z Z 00 m,m1 (cid:20) m′,m′1 (cid:21)
mX1,m′1 0 0
Integrationovervariables θ,φ in the aboveformula is performedwith the use offormulae (4.6.2), (4.1.12)and (4.2.7)
from the reference [6]. They lead to expression
πdθ 2πdφsinθD(l1)(0,θ,φ)D(l) (0,θ,φ) D(l′) (0,θ,φ) ∗
Z Z 00 m,m1 (cid:20) m′,m′1 (cid:21)
0 0
=4π( 1)m′−m′1 l l′ l1 l l′ l1 ,
− (cid:18)m −m′ 0 (cid:19)(cid:18)m1 −m′1 0 (cid:19)
5
which contains Wigner 3-j symbols. Taking into account the above integral in expression (9) yields
Gˆlm,l′m′(kˆez)=4πl1|=Xl+|ll−′|l′|mX1=l −lmX′1=l′−l′(−1)m′−m′1(2l1+1)(cid:18)ml −lm′ ′ l01 (cid:19)(cid:18)ml1 −lm′ ′1 l01 (cid:19)×
∞
(−i)l1Z dRR2jl1(kR)Glm1,l′m′1(Rˆez). (10)
0
In this way the three-dimensional Fourier transform is reduced to the one dimensional Hankel transform [9].
ItisconvenienttointroducedifferentrepresentationofmatrixG. Namely,insteadofGlm,l′m′(Rˆez)wecanconsider
gj (R) defined in the following way
l,l′
l l′ j
glj,l′(R)=(2j+1) (−1)m(cid:18)m m′ 0(cid:19)Glm,l′m′(Rˆez). (11)
mX,m′ −
with the inverse transformation
l+l′ l l′ j
Glm,l′m′(Rˆez)=δm,m′(−1)m (cid:18)m m′ 0 (cid:19)glj,l′(R). (12)
j=X|l−l′| −
Let us notice that only integer j which satisfy condition l l′ j l +l′ needs to be considered in the above
equations. It follows from properties of Wigner 3-j symbol|s.−Th|e≤mult≤ipole matrix elements Glm,l′m′(R) for any R
are then related to gj (R) by equation
l,l′
Glm,l′m′(R(R,θ,φ))=Xm1 (−1)m′+m1j=Xl|+l−l′l′|(cid:18)ml −lm′ ′ mj1 (cid:19)D−(jm) 1,0(0,−θ,−φ)glj,l′(R), (13)
which follows from equation (8) and relation (12).
Different representations of the multipole matrix Glm,l′m′(Rˆez) in positional space given by equations (11) and
(12) can be similarly introduced also in Fourier space, i.e.
l l′ j
g˜lj,l′(k)=(2j+1) (−1)m(cid:18)m m′ 0(cid:19)Gˆlm,l′m′(kˆez). (14)
mX,m′ −
l+l′ l l′ j
Gˆlm,l′m′(kˆez)=δm,m′(−1)m (cid:18)m m′ 0 (cid:19)g˜lj,l′(k). (15)
j=X|l−l′| −
By equations (10), (14), and (12) and orthogonality of Wigner 3-j symbols [6], we get that in the new basis the
Fourier transform of the multipole matrix is expressed in the following form
∞
g˜j (k)=4π( i)j dRR2j (kR)gj (R). (16)
l,l′ − Z j l,l′
0
Inthiswaythethree-dimensionalFouriertransformgivenbyexpression(2)ofamultipolematrixsatisfyingrotational
symmetry (5) is reduced to relevant Hankel transform given by expression (16).
Calculations performed in this section can be repeated with very minor modification in order to reduce the inverse
three-dimensional Fourier transform of a multipole matrix to the one-dimensional Hankel transform. We omit the
derivationgivingonlytheresultfortheinverseFouriertransformofamultipolematrixsatisfyingrotationalsymmetry,
1 ∞
gj (R)= ij dk k2j (kR)g˜j (k). (17)
l,l′ 2π2 Z j l,l′
0
6
V. MULTIPOLE MATRIX IN POSITIONAL SPACE
To calculate the inverse Fourier transform of Gˆ(k) with the use of equation (17), g˜j (k) is needed. We
h ilm,l′m′ ll′
calculate it by means of transformation (14) and expression (3) for Gˆlml′m′(kˆez). The calculations yields
g˜j (k)=4π( i)−l+l′(2j+1)[(2l+1)(2l′+1)]1/2al+l′+2 l l′ j jl(ka)jl′(ka). (18)
l,l′ − (cid:18)0 0 0(cid:19) k2
It is worth noting that the Wigner 3 j symbol is not vanishing only when the l,l′,j satisfy triangular inequality,
l l′ j l+l′ . In calculations of−the above formula we used the following property of the spherical harmonics,
| − |≤ ≤| |
Y (ˆe )=δ ((2j+1)/4π)1/2.
jm z m,0
With the above expression for g˜j (k), equation (17) yields
l,l′
∞
glj,l′(R)=µj,l,l′al+l′+2Z dk jj(kR)jl(ka)jl′(ka) (19)
0
with
µj,l,l′ = 2(−i)−l+πl′+j(−1)j (2j+1)[(2l+1)(2l′+1)]1/2(cid:18)0l l0′ 0j (cid:19). (20)
ItdemandstocalculateintegralofthreesphericalBesselfunctions. Thathasalreadybeenconsideredintheliterature
[10]. For R 2a, the integral is given as follows
≥
∞
Z dk jl1(kR)jl(ka)jl′(ka)
0
π3/2 a l+l′+1 Γ 1 +l+l′
= 8a δl+l′,l1(cid:16)R(cid:17) Γ 3(cid:0)+2l Γ 3 +(cid:1)l′ , (21)
2 2
(cid:0) (cid:1) (cid:0) (cid:1)
with Euler Gamma function Γ(x). Whereas for R 2a we have
≤
∞
dk jj(kR)jl(ka)jl′(ka)
Z
0
π3/2R l l′ 1+l l′ 1 l+l′ 2+l+l′ 1 3 j 4+j R2
= 2a aαj,l,l′ 4F3(cid:18)− 2− , 2− , −2 , 2 ;2, −2 , 2 ;4a2(cid:19)
π3/2 R j j l l′ 1 j+l l′ j+l′ l l+l′+j+1 1+j j 3 R2
+ 2a (cid:18)a(cid:19) βj,l,l′ 4F3(cid:18) − −2 − , 2− , 2 − , 2 ; 2 ,2,2 +j;4a2(cid:19)
π3/2R2 l′+l+1 l l′ l′ l l+l′+1 3 j j+1 R2
− 2a a2γj,l,l′ 4F3(cid:18)1− 2 ,1+ −2 ,1+ 2− ,1+ 2 ;2,2− 2,2+ 2 ;4a2(cid:19) (22)
with coefficients αj,l,l′, βj,l,l′, and γj,l,l′ given by
Γ j−1
αj,l,l′ =2−5/2Γ 1+l′−l Γ(cid:0)1+2l−(cid:1)l′ Γ j+4 ,
2 2 2
(cid:0) (cid:1) (cid:0) Γ((cid:1)1 (cid:0)j)Γ(cid:1) 1+l+l′+j
βj,l,l′ =2−3/2Γ 1+ 1+l+l′−j Γ 1−+ l′−l−(cid:16)j Γ2 1+(cid:17)l−l′−j Γ 3 +j ,
2 2 2 2
(cid:16) (cid:17) (cid:16) (cid:17) (cid:16) (cid:17) (cid:0) (cid:1)
Γ j 1
γj,l,l′ =−2−7/2Γ l′−l Γ (cid:0)l−2l′−Γ(cid:1)2+ j+1 . (23)
2 2 2
(cid:0) (cid:1) (cid:0) (cid:1) (cid:0) (cid:1)
Symbol F stands for hypergeometric function
4 3
∞ (α ) (α ) (α ) (α ) xk
F (α ,α ,α ,α ;β ,β ,β ;x)= 1 k 2 k 3 k 4 k
4 3 1 2 3 4 1 2 3
(β ) (β ) (β ) k!
Xk=0 1 k 2 k 3 k
7
where (α) denotes Pochhammer symbol defined by
k
(α) =1 for k =0
0 .
(α) =α(α+1)...(α+k 1) for k =1,2,...
k −
For nonoverlapping configurations R > 2a, the multipole matrix elements are given by formulae (19), (20), and
(21). Therefore in this regime gj (R) is proportional to 1/Rl+l′+1.
ll′
For overlapping configurations, R < 2a, the multipole matrix elements are given by formulae (19), (20), and (22)
with equations (23). For small l and l′, analysis of the poles in Euler gamma functions and in Pochhammer symbols
(which may also be expressed by Euler gamma functions) reveals that for overlapping configurations gj (R) is a
ll′
polynomial with respect to R. We observe that the degree of the polynomial is l+l′+1. For example
( 2a+R)2(4a+R)
g0 (R<2a)= − , (24)
1,1 − 16√3
4a2R+R3 2
g3 (R<2a)= 7 − . (25)
2,3 − (cid:0) 256√3 (cid:1)
VI. SUMMARY
In this article we calculate the multipole matrix elements of Green function for Laplace equation. The elements
are defined by equation (1). The expression for nonoverlapping configurations, i.e. for R > 2a, is known in the
literatureandcanbe inferrede.g. fromthe reference[7]. The new contributionis the calculationofexpression(1)for
overlapping configurations, i.e. for R < 2a. In this case one can find related considerations for the lowest multipole
numbers [8].
It is worth mentioning that the method which we use to calculate multipole matrix for Laplace equation can be
generalizedforothercases,becausetheconsiderationsrelyonsimplificationoftheFouriertransform. Forexample,the
methodcanbegeneralizedtothe multipolematrixelementsofGreenfunctionforStokesequationsinhydrodynamics
[11–13]. In this case overlapping configurations of multipole Green function for the lowest multipoles have been
considered recently in the literature [14].
ACKNOWLEDGMENTS
K.M. has been supported by MNiSW grant IP2012 041572,and, at the earlier stage of the research, also acknowl-
edged support by the Foundation for Polish Science (FNP) through the TEAM/2010-6/2project, co-financed by the
EU European Regional Development Fund.
[1] M. Sahimi, Heterogeneous materials (Springer, 2003) ISBN 0387001662.
[2] D.Griffiths and R. College, Introduction to electrodynamics, Vol.3 (prentice Hall New Jersey, 1999).
[3] S.KimandS.Karrila,Microhydrodynamics: PrinciplesandSelectedApplications (Butterworth-HeinemannBoston,1991).
[4] M. H.Sadd, Elasticity: theory, applications, and numerics (Academic Press, 2014).
[5] K.Hinsen and B. Felderhof, Physical Review B 46, 12955 (1992).
[6] A.Edmonds, Angular momentum in quantum mechanics (Princeton Univ Pr, 1996).
[7] S.-Y.Sheu,S.Kumar, and R. Cukier, Physical Review B 42, 1431 (1990).
[8] B. Felderhof, G. Ford,and E. Cohen, Journal of Statistical Physics 33, 241 (1983), ISSN0022-4715.
[9] N.Bleistein and R. Handelsman, Asymptotic expansions of integrals (Harcourt College Pub, 1975) ISBN 0030835968.
[10] A.Prudnikov,Y.A. Bryčkov,and O. Maričev(1983).
[11] M.L.Ekiel-JeżewskaandE.Wajnryb,“Precisemultipolemethodforcalculatinghydrodynamicinteractionsbetweenspher-
icalparticlesinthestokesflow,in: Theoreticalmethodsformicroscaleviscousflows,transworldresearchnetwork,” (2009).
[12] B. Cichocki, R.Jones, R.Kutteh,and E. Wajnryb,The Journal of Chemical Physics 112, 2548 (2000).
[13] B. Felderhof and R. Jones, Journal of Mathematical Physics 30, 339 (1989).
[14] E. Wajnryb,K. A.Mizerski, P. J. Zuk, and P. Szymczak,Journal of Fluid Mechanics 731, R3(2013).