A Series from an integral representation

# Exact Fourier expansion in cylindrical coordinates for the three-dimensional Helmholtz Green function

## Abstract

A new method is presented for Fourier decomposition of the Helmholtz Green Function in cylindrical coordinates, which is equivalent to obtaining the solution of the Helmholtz equation for a general ring source. The Fourier coefficients of the Helmholtz Green function are split into their half advanced+half retarded and half advancedhalf retarded components. Closed form solutions are given for these components in terms of a Horn function and a Kampé de Fériet function, respectively. The systems of partial differential equations associated with these two-dimensional hypergeometric functions are used to construct a fourth-order ordinary differential equation which both components satisfy. A second fourth-order ordinary differential equation for the general Fourier coefficent is derived from an integral representation of the coefficient, and both differential equations are shown to be equivalent. Series solutions for the various Fourier coefficients are also given, mostly in terms of Legendre functions and Bessel/Hankel functions. These are derived from the closed form hypergeometric solutions or an integral representation, or both. Numerical calculations comparing different methods of calculating the Fourier coefficients are presented.

## 1 Introduction and overview

The inhomogeneous Helmholtz wave equation is

 (∇2+β2)Φ(β,r)=ρ(r) (1)

and this has the well known free-space retarded Green function [1, p. 284]

 GH(β,r−r′)=−exp(iβ∣r−r′∣)4π∣r−r′∣ (2)

where is a field point, is a source point and is the wave number, here considered to be a general complex number. The free-space Green function (2) is restricted to values of such that as . For general dispersive waves with where and are real, then is a condition for this to hold. In the limit as these equations reduce to the Poisson equation and its corresponding Green function . The general retarded solution of the Helmholtz equation at a field point for a general source density , subject to the boundary condition that as , is given in terms of the Green function as

 Φ(β,r)=\diiintGH(β,r−r′)ρ(r′)d3r′ (3)

where the volume integral is to be taken over all regions of space where the source density is non zero.

Many problems of practical interest have some element of axial symmetry and are best treated in cylindrical coordinates , the Cartesian components of being related to the cylindrical components by . It follows immediately from this relation that the distance between a source point and a field point is given by

 ∣r−r′∣=√r2+r′2+(z−z′)2−2rr′cos(ϕ−ϕ′). (4)

The solution for when is a general circular ring source is of particular interest, with applications such as circular loop antennas [2], [3], [4], the acoustics of rotating machinery [5] and acoustic and electromagnetic scattering [6]. For the simpler Poisson equation most of the analytical solutions found in the literature for cylindrical geometry are either ring source solutions or can be easily constructed from them by integration or summation. Examples are gravitating rings and disks, ring vortices and vortex disks, and circular current loops and solenoids.

The source density for a thin circular ring of radius located in the plane  is of the form

 ρc(r,R,Z)=f(ϕ)δ(r−R)δ(z−Z) (5)

where is the angular distribution of the source strength around the ring. This can be most conveniently described by a Fourier series of the form

 \limits is allowed only on operators (6)

where the Fourier coefficients and are given by [7, p. 1066]

 \limits is allowed only on operators (7)
 \limits is allowed only on operators (8)

From equation (4), the Green function (3) is even in the variable , where is the angular coordinate of and is the angular coordinate of . It is convenient to exploit this symmetry when substituting equations (5) and (6) into (3). From the identity we obtain

 \limits is allowed only on operators
 \limits is allowed only on operators (9)

and on substituting equations (5) and (9) into (3) and performing the volume integration, the odd terms proportional to in equation (9) do not contribute to the solution as is even in . The remaining integrals from the even terms can be calculated over the reduced interval from to . This gives the solution of the Helmholtz equation for a circular ring source with general in the form

 Φc(β,r,R,Z)=−a02G0H(β,r,R,z−Z)
 \limits is allowed only on operators (10)

where

 \limits is allowed only on operators (11)

and where the explicit dependence of the solution on the constant ring parameters and has been introduced in these definitions. Introducing the Neumann factor such that for and for , and defining allows (10) to be expressed more concisely as

 \limits is allowed only on operators (12)

Apart from a constant factor, the terms in (12) are also the coefficients in the Fourier expansion of the Green function (2) itself, when the source point is given by . From equations (6), (7) and (8) this is given by

 \limits is allowed only on operators (13)

Thus the solution of the Helmholtz equation for a general ring source can be constructed directly from the coefficients in the Fourier expansion of the Green function (3). This provides in large measure the motivation to analytically construct the Fourier series for the Helmholtz Green function.

For the Poisson equation with the corresponding Fourier expansion of the Green function has already been given in closed form as [8]:

 \limits is allowed only on operators (14)

where

 ω=r2+R2+(z−Z)22rR (15)

is a toroidal variable such that and the are the Legendre functions of the second kind and half-integral degree, which are also toroidal harmonics. The Fourier expansion given by equations (14) and (15) can be obtained immediately by writing the Green function (2) for in the form

 GP(r−r′)=−14π√rR√2ω−2cos(ϕ−ϕ′) (16)

where is given by (15), and noting that the function has the simple integral representation [7, eqn 8.713]

 \limits is allowed only on operators (17)

An alternative derivation of (14) employs the Lipschitz integral [7, eqn 6.611 1]

 \limits is allowed only on operators (18)

and Neumann’s addition theorem [9, eqn11.2 1]

 \limits is allowed only on operators (19)

to obtain the well known eigenfunction expansion

 GP(r−r′)=−14π×
 \limits is allowed only on operators (20)

This reduces to (14) on employing the integral [7, eqn 6.612 3],[9, eqn 13.22]:

 \limits is allowed only on operators (21)

The generalization of (20) for the Helmholtz case is also well known [10, p. 888]

 GH(β,r−r′)=−i4π×
 \limits is allowed only on operators (22)

This can be similarly obtained from Neumann’s theorem by employing the integral [7, eqn 6.616 2]

 \limits is allowed only on operators (23)

instead of the Lipschitz integral. Equation (22) gives the Fourier coefficients of the Helmholtz Green function in the form

 \limits is allowed only on operators (24)

This reduces to (20) in the limit as but unfortunately the integral in (24) is not given in standard tables for . Numerical evaluation of this integral requires care, as the integrand is oscillatory and singular in an infinite range of integration, though the integrand tends exponentially to zero as . Equation (11) is a convenient alternative numerical evaluation of the Fourier coefficients, provided is not too large.

The integrals (11) and (24) contain the additional parameter which is not contained in (17) and (21). As a consequence of this, the closed form generalization of (14) for the Helmholtz case involves two-multidimensional Gaussian hypergeometric series, and the main purpose of this article is to present these solutions and various related results. The core idea leading to the solution is expansion of the exponential in (2) as the absolutely convergent power series [4]

 \limits is allowed only on operators (25)

where

 Missing or unrecognized delimiter for \right (26)

Hence

 \limits is allowed only on operators (27)

where

 \limits is allowed only on operators (28)

The integral (28) can be evaluated as a series by binomial expansion and this gives a double series for the Fourier coefficient . The expansion of (28) gives an infinite number of terms for even and a finite number of terms for odd. These two cases are best treated separately and it is therefore convenient to split the summation over in (25) into odd and even terms. This is equivalent to splitting the Green function (2) such that

 GH(β,r−r′)=Λ+(β,r−r′)+Λ−(β,r−r′) (29)

where

 Λ+(β,r−r′)=−18π(exp(iβ∣r−r′∣)∣r−r′∣+exp(−iβ∣r−r′∣)∣r−r′∣) (30)

is the half advanced+half retarded Green function and

 Λ−(β,r−r′)=−18π(exp(iβ∣r−r′∣)∣r−r′∣−exp(−iβ∣r−r′∣)∣r−r′∣) (31)

is the half advancedhalf retarded Green function. The corresponding Fourier coefficients are split in the same manner such that

 Λm+(β,r,R,z−Z)=12(GmH(β,r,R,z−Z)+GmH(−β,r,R,z−Z)) (32)
 Λm−(β,r,R,z−Z)=12(GmH(β,r,R,z−Z)+GmH(−β,r,R,z−Z)). (33)

For real , splitting the Green function in this way is equivalent to dividing it into its real and imaginary parts, but this is not the case for general complex . It is shown in Section 2 that the Fourier coefficients in (32) and (33) are given respectively by

 ¯Λm+(rR,γ,k)=
 Γ(m+1/2)m!√πrR(k2)2m+1H3(m+1/2,m+1/2,2m+1,k2,γ24) (34)

and

 ¯Λm−(rR,γ,k)=i√rR(2m+1)!(γk2)2m+1×
 F0:1;01:1;0[−:m+3/2:m+1/2;2m+1;−;−;γ2k24,−γ24] (35)

where

 k=√4rR(r+R)2+(z−Z)2 (36)
 γ=β√(r+R)2+(z−Z)2 (37)

and

 ¯Λm±(rR,γ,k)≡Λm±(β,r,R,z−Z). (38)

The variable is the usual modulus contained in elliptic integral solutions of elementary ring problems and is related to the toroidal variable by

 ω=2−k2k2. (39)

The function in equation (34) is one of the standard Horn functions [11, eqn 5.7.1 31] and is equivalent to the double hypergeometric series

 H3(m+1/2,m+1/2,2m+1,k2,γ24)=
 \limits is allowed only on operators (40)

The Kampé de Fériet function [12, p. 27] in (35) is equivalent to the double hypergeometric series

 F0:1;01:1;0[−:m+3/2:m+1/2;2m+1;−;−;γ2k24,−γ24]=
 \limits is allowed only on operators (41)

The integral (28) can also be evaluated using an integral representation for the associated Legendre function of the first kind, and it is shown in Appendix A that this gives the series expansion:

 Missing or unrecognized delimiter for \left
 \limits is allowed only on operators (42)

where

 λ=β√2rR (43)

and

 ^GmH(rR,λ,ω)≡GmH(β,r,R,z−Z). (44)

The Legendre function in equation (42) reduces to an associated Legendre polynomial for odd . The series in (42) can be split into even and odd terms such that

 ^GmH(rR,λ,ω)=12(^Λm+(rR,λ,k)+^Λm−(rR,λ,ω)) (45)

where

 ^Λm±(rR,λ,ω)≡Λm±(β,r,R,z−Z) (46)

and it is shown in Appendix A that the even and odd series can be expressed respectively as:

 \limits is allowed only on operators (47)
 \limits is allowed only on operators (48)

In equation (48) the Legendre function is purely imaginary for real . In the static limit as then and from the gamma function identity [7, eqn 8.334 2]

 Γ(1/2−m)Γ(1/2+m)=(−1)mπ % [m∈N0] (49)

then equation (47) reduces to

 ^Λm+(rR,0,ω)=Qm−1/2(ω)π√rR [m∈N0] (50)

as it must do for consistency with (14).

The solutions in terms of two-dimensional hypergeometric functions defined by Equations (34)-(38) and (40)-(41) can be summed over either index to give the solutions as series of special functions. It is shown in Section 3 that summation over the index in equation (40) gives equation (47), exactly as given by the integral representation. However, summation over the index in equation (41) gives instead the series solution

 ¯Λm−(rR,γ,k)=i√πm!√rR(γk4)2m+1×
 \limits is allowed only on operators (51)

A hypergeometric identity to reduce the hypergeometric function in equation (51) to other well-known special functions does not seem to be available in standard tabulations. It might nevertheless be conjectured that (51) could somehow be reducible to equation (48), but this is not in fact the case. It is easily verified numerically that although equations (48) and (51) both converge rapidly to the same limit, the individual terms do not match. Hence, equation (51) is a distinct series from equation (48). It is also shown in Section 3 that summation over the index in equations (40) and (41) gives the Bessel function series:

 \limits is allowed only on operators (52)
 ¯Λm−(rR,γ,k)=i2√rR∞∑n=0Γ(n+m+1/2)Γ(n+2m+1)n!(γk22)n+m+1/2Jn+m+12(γ) (53)

and these two series can be conveniently combined to give a series of Hankel functions of the first kind:

 \limits is allowed only on operators (54)

From the solutions (34) and (35) it can be seen that dimensionless Fourier coefficients defined by

 gm±(k2,γ2/4)≡√rR¯Λm±(rR,γ,k) (55)

depend only on the two dimensionless variables and . The functions are given explicitly by equations (34) and (35) as:

 gm+(x,y)=Γ(m+1/2)22m+1m!√πxαH3(α,α,2α,x,y) (56)

and

 gm−(x,y)=i(2m+1)!(xy)αF0:1;01:1;0[−:α+1:α;2α;−;−;xy,−y] (57)

where

 α=m+1/2. (58)

Two-dimensional hypergeometric series such as (56) and (57) are associated with pairs of partial differential equations [11, section 5.9] and these can be used to construct ordinary differential equations for with fixed and as the independent variable. It is shown in Section 4 that for constant the coefficients both satisfy the same fourth-order ordinary differential equation in :

 (1−x)x4d4gm±dx4+(6−9x)x3d3gm±dx3+(α(1−α)+6−18x−xy+2y)x2d2gm±dx2+
 (2α(1−α)−2x(3−y))xdgn±dx+(y2+2α(1−α)y−3α(α+1)(α+2))gm±=0. (59)

In Section 5 an integral representation is derived for

 ^ym(λ,ω)≡π√2rR^GmH(rR,λ,ω) (60)

and this is used to derive a fourth-order ordinary differential equation for in terms of :

 (1−ω2)d4^ymdω4−6ωd3^ymdω3+(m2−λ2ω2−254)d2^ymdω2−λ2d^ymdω−(λ24)2^ym=0. (61)

In the static limit as , equation (61) reduces to:

 d2dω2[(1−ω2)d2^ymdω2−2ωd^ymdω+(m2−14)^ym]=0 (62)

where

 (1−ω2)d2^ymdω2−2ωd^ymdω+(m2−14)^ym=0 (63)

is Legendre’s equation of degree [7, eqn 8.820]. It is also shown in Section 5 that the differential equations (59) and (61), obtained by quite different routes, are equivalent. The special functions used in the analysis are given in Table 1.

Recurrence relations for the Fourier coefficients for the Helmholtz equation were investigated by Matviyenko [6], but the closed form solutions and differential equations presented here appear to be new. Werner [3] presented an expansion of the Fourier coefficient as a series of spherical Hankel functions, superficially similar to equation (54), but the two expansions are distinct. The two-dimensional hypergeometric series approach applied here to obtain the Fourier expansion for the Helmholtz Green function has recently been applied to obtain the Fourier expansion in terms of the amplitude for the Legendre incomplete elliptic integral of the third kind [13].

The numerical performance of the various expressions for the Fourier coefficients was investigated using Mathematica [14] and this is examined in Appendix C.

## 2 Solution in terms of two-dimensional hypergeometric series

The power series expansion (27) for the Fourier coefficient can be expressed in the form

 \limits is allowed only on operators (64)

where

 GmH(β,r,R,z−Z)=¯GmH(rR,γ(β,r,R,z−Z),k(r,R,z−Z)) (65)

and and are defined by (36) and (37). The term in (64) can be expanded binomially to give as a double series containing integrals of the form

 \limits is allowed only on operators (66)

This integral is given by Gradshteyn and Ryzhik [7, eqns 3.631 8,12] in a form which can be recast as

 Im,p=(−1)mπ22p+1(2p+1)1B(p+m+1,p−m+1) for p⩾m (67)
 Im,p=0 for p

and expressing the beta function in (67) in terms of gamma functions and employing the duplication theorem [7, eqn 8.335 1]

 Γ(2x)=22x−1√πΓ(x)Γ(x+1/2) (69)

gives after some reduction the alternative form

 Im,p=(−1)m√π2Γ(p+1/2)p!Γ(p+m+1)Γ(p−m+1) for p⩾m. (70)

The binomial expansion of (64) gives an infinite series for zero or even and a finite sum for odd, and these two cases must be treated separately. It is therefore convenient to split the series for such that

 ¯GmH(rR,γ,k)=¯Λm+(rR,γ,k)+¯Λm−(rR,γ,k) (71)

and on employing equation (69) the divided series are given by

 \limits is allowed only on operators (72)
 \limits is allowed only on operators (73)

Binomial expansion of the integrals in equations (72) and (73) gives respectively

 \limits is allowed only on operators (74)
 \limits is allowed only on operators (75)

and employing the explicit formula (70) for in (74) and (75) gives respectively

 ¯Λm+(rR,γ,k)=k2π√rR×
 \limits is allowed only on operators (76)
 ¯Λm−(rR,γ,k)=ikγ4√rR×
 \limits is allowed only on operators (77)

where the gamma identity (49) has been used to simplify equation (76).

### 2.1 The series for ¯Λm+(rR,γ,k)

The substitution in equation (76) yields after some reduction the double series:

 ¯Λm+(rR,γ,k)=Γ(m+1/2)m!√πrR(k2)2m+1×
 \limits is allowed only on operators (78)

where in (78)

 (a)p≡Γ(a+p)Γ(a) (79)

is the Pochhammer symbol. The double hypergeometric function in (78) can be identified as one of the confluent Horn functions [11, eqn 5.7.1 31] and hence

 ¯Λm+(rR,γ,k)=Γ(m+1/2)m!√πrR(k2)2m+1×
 H3(m+1/2,m+1/2,2m+1,k2,γ24). (80)

The convergence condition given in [11, eqn 5.7.1 31] for the double series in equation (78) is , which always holds. The order of summation in equation (78) can be reversed, but the order of the arguments and in equation (80) cannot be exchanged.

### 2.2 The series for ¯Λm−(rR,γ,k)

Equation (77) can be converted to a doubly infinite series by reversing the order of summation, which gives

 ¯Λm−(rR,γ,k)=ikγ4√rR×
 \limits is allowed only on operators (81)

The substitution in (81) gives

 ¯Λm−(rR,γ,k)=ikγ4√rR×
 \limits is allowed only on operators (82)

and the further substitution in (82) gives

 ¯Λm−(rR,γ,k)=i2√rR×(γ2k24)m+1/2×
 \limits is allowed only on operators (83)

Expressing equation (83) in terms of Pochhammer symbols gives after some reduction the double hypergeometric series

 ¯Λm−(rR,γ,k)=i√rR(2m+1)!(γk2)2m+1×
 \limits is allowed only on operators (84)

and this can be expressed as a Kampé de Fériet function as defined by Srivastava and Karlsson [12, p. 27]:

 Fp:q;rl:m;n[{ap}:{αl}:{bq};{βm};{cr};{γn};x,y]=
 \limits is allowed only on operators (85)

In the definition (85), and and so on, are the lists of the arguments of the Pochhammer symbols of the various types which appear in the products on the right-hand side of the equation. If a list has no members, it is represented by a hyphen. Comparing (84) with (85) gives

 ¯Λm−(rR,γ,k)=i√rR(2m+1)!(γ2k24)m+1/2×
 F0:1;01:1;0[−:m+3/2:m+1/2;2m+1;−;−;γ2k24,−γ24]. (86)

## 3 Consequences of the hypergeometric formulas

In the static limit as then equation (34) reduces to

 ¯Λm+(rR,0,k)=Γ(m+1/2)m!√πrR(k2)2m+12F1(m+1/2,m+1/2;2m+1;k2) (87)

and from (37) and the standard hypergeometric identity [15, eqn 7.3.1 71]:

 2F1(a,b;2b;z)=22b√πΓ(b+1/2)Γ(2b−a)z−b(1−z)(b−a)/2exp(iπ(a−b))Qb−ab−1(2z−1) (88)

this reduces to

 ¯Λm+(rR,0,k)=1π√rRQm−1/2(2−k2k2) (89)

in agreement with equations (14) and (37).

### 3.1 Fourier coefficients as series of special functions

The double series given by equation (34) can be summed with respect to either the index or the index in the definition (40). Summing with respect to in (40) gives a series of Bessel functions of the second kind:

 \limits is allowed only on operators (90)

where the gamma function identity (49) and the Bessel function identity

 Yν(γ)=1sinνπ[cos(νπ)−J−ν(γ)] (91)

have been employed to obtain equation (90). Summing instead over the index in (40) gives the alternative series

 ¯Λm+(rR,γ,k)=1m!√πrR(k2)2m+1×
 \limits is allowed only on operators (92)

This can be reduced using (88) and (49) to give:

 ¯Λm+(rR,γ,k)=(−1)m√rR×
 \limits is allowed only on operators (93)

The dimensionless variables and defined by equations (37) and (43) respectively are related by

 γ=√2λk (94)

and substituting this equation and equation (39) in equation (93) gives immediately equation (47).

Summing with respect to in equation (83) gives the Bessel series

 ¯Λm−(rR,γ,k)=i2√rR∞∑p=0Γ(p+m+1/2)Γ(p+2m+1)p!(γk22)p+m+1/2Jm+12+p(γ) (95)

The corresponding summation over the index gives

 ¯Λm−(rR,γ,k)=i√πm!√rR(γk4)2m+1×
 \limits is allowed only on operators (96)

There seems to be no hypergeometric transformation listed in standard tables suitable for directly reducing the hypergeometric function in this equation.

Equations (90) and (95) can be conveniently combined to give a series of Hankel functions of the first kind:

 \limits is allowed only on operators (97)

where .

## 4 Differential equations from the two-dimensional hypergeometric solutions

Erdélyi et. al. [11, section 5.9] tabulate the partial differential equations satisfied by all the functions in Horn’s list. They employ the notation reproduced below for the various partial derivatives

 p=∂z∂x; q=∂z∂y;% r=∂2z∂x2; s=∂2z∂x∂y; t=∂2z∂y2

where is any function on the list. Each function in the list satisfies two partial differential equations, and unfortunately those given in [11, 5.9 34] for contain typographical errors. The correct equations can be shown by the methods given in [11, section 5.7] to be:

 x(1−x)r+xys+[δ−(α+β+1)x]p+βyq−αβz=0 (98)
 yt−xs+(1−α)q+z=0 (99)

and for the particular case considered here we have and so that (98) reduces to

 x(1−x)r+xys+[2α−(2α+1)x]p+αyq−α2z=0. (100)

Similar equations can also be derived for the Kampé de Fériet function defined by equation (41). For the definition

 ¯z(α,u,v)≡F0:1;01:1;0[−:α+1:α;2α;−;−;u,v] (101)

and the notation

 ¯p=∂¯z∂u; ¯q=∂¯z∂v; ¯r=∂2¯z∂u2; ¯s=∂2¯z∂u∂v; ¯t=∂2¯z∂v2

then the equations corresponding to (99) and (100) can be shown to be

 ¯z=v¯t+u¯s+(α+1)¯q (102)
 ¯z=2α¯p+u¯r+¯q+v¯t. (103)

### 4.1 Fourth-order differential equation from the Horn Function H3(α,β,δ,x,y)

Writing

 z(α,x,y)=H3(α,α,2α,x,y) (104)

then (56) can be expressed as

 z(α,x,y)=22m+1m!√πΓ(m+1/2)x−αgm+(x,y). (105)

The ordinary differential equation for in terms of can be derived by first obtaining the corresponding differential equation for from (99) and (100) and then substituting (105) into this equation. Although straightforward in principle, this procedure is rather intricate in practice, and only the essential elements of the derivation are given below.

#### Differential equation for z(α,x,y)

It is convenient to define a differential operator such that:

 D≡xddx (106)

which has by definition the properties:

 Dx=x (107)
 Dy=−y (108)
 D(xy)=0 (109)
 Dz=xp−yq (110)
 Dp=xr−ys (111)
 Dq=xs−yt. (112)

Combining equations (99) and (112) gives the equation

 Dq=z+(1−α)q (113)

and applying the operators and to this equation gives respectively:

 Ds=p−αs (114)
 Dt=q+(2−α)t. (115)

Eliminating the variable gives a system of four coupled equations:

 Dz+yq=xp (116)
 (1−x)xDp=−xys+[(2α+1)x−2α]xp−αxyq+α2xz (117)
 Dq=z+(1−α)q (118)
 xDs=xp−αxs. (119)

Eliminating gives a system of 3 equations:

 D2z−xD2z+(2α−1)Dz−2αxDz+yz−xyz−α2xz=−xys+(1−α)yq (120)
 Dq=z+(1−α)q (121)
 xDs=Dz+yq−αxs. (122)

Eliminating gives the two equations:

 D3z−xD3z−(3α+1)xD2z+(3α−1)D2z−α(3α+2)xDz+2yDz
 +α(2α−1)Dz−xyDz−α2(α+1)xz+2(α−1)yz−αxyz=−y2q (123)
 Dq=z+(1−α)q. (124)

Eliminating gives finally the fourth-order equation:

 (