1. Introduction
Internal gravity waves in density-stratified fluids, inertial waves in rotating fluids and inertia–gravity waves in rotating stratified fluids, all share a common pattern under localized monochromatic excitation: a St Andrew's cross in two dimensions, and a double cone in three dimensions. The first experimental studies of this pattern used oscillating bodies to generate the waves, such as a circular cylinder (Mowbray & Rarity Reference Mowbray and Rarity1967) and a sphere (McLaren et al. Reference McLaren, Pierce, Fohl and Murphy1973). The first theoretical studies were modelled after these experiments, and considered cylinders either circular (Hurley Reference Hurley1972; Appleby & Crighton Reference Appleby and Crighton1986) or elliptic (Hurley Reference Hurley1997; Hurley & Hood Reference Hurley and Hood2001), and spheres (Hendershott Reference Hendershott1969; Appleby & Crighton Reference Appleby and Crighton1987; Voisin Reference Voisin1991; Rieutord, Georgeot & Valdetarro Reference Rieutord, Georgeot and Valdetarro2001; Davis Reference Davis2012) or spheroids (Krishna & Sarma Reference Krishna and Sarma1969; Sarma & Krishna Reference Sarma and Krishna1972; Lai & Lee Reference Lai and Lee1981). The free oscillations of a sphere, displaced from its neutral buoyancy level then released, were also considered (Larsen Reference Larsen1969).
The theory combined coordinate stretching and analytic continuation; an equivalent formulation involving the Laplace transform for impulsively started oscillations was also used. The problem was solved first in the frequency range where the waves are evanescent and their equation elliptic, by stretching the coordinates so as to turn the wave equation into the Laplace equation when the Boussinesq approximation is made, and a Helmholtz equation otherwise; the solution was then continued analytically to the frequency range where the waves are propagating and their equation hyperbolic, based on the causality requirement that, for time dependence as $\exp (-\mathrm {i}\omega t)$ say, with t the time and
$\omega$ the frequency, the solution be analytic in the upper half of the complex
$\omega$-plane.
Several origins have been invoked for this approach: Hendershott (Reference Hendershott1969) traced it back to Stewartson (Reference Stewartson1952), for the formation of Taylor columns; Appleby & Crighton (Reference Appleby and Crighton1986) traced it back to Hurley (Reference Hurley1972), who combined it with conformal mapping to devise a general approach of two-dimensional monochromatic internal waves; and Voisin (Reference Voisin1991) traced it back to Pierce (Reference Pierce1963), for the Green's function of acoustic–gravity waves. But Rieutord et al. (Reference Rieutord, Georgeot and Valdetarro2001) pointed out that the approach is actually much older, having been introduced first by Bryan (Reference Bryan1889) for the determination of the inertial modes of a rotating spheroid, a problem to which it is still applied to this day (Rieutord & Noui Reference Rieutord and Noui1999; Ivers, Jackson & Winch Reference Ivers, Jackson and Winch2015; Backus & Rieutord Reference Backus and Rieutord2017).
As pointed out by Ermanyuk (Reference Ermanyuk2002), the approach also allows the determination of the added mass of a rigid body oscillating in a stratified fluid, based on the added mass of the stretched version of the same body oscillating in a homogeneous fluid. The outcome has been compared with experiment for a spheroid (Ermanyuk Reference Ermanyuk2002), a sphere (Ermanyuk & Gavrilov Reference Ermanyuk and Gavrilov2003), a cylinder with diamond-shaped (Ermanyuk & Gavrilov Reference Ermanyuk and Gavrilov2002) or circular (Brouzet et al. Reference Brouzet, Ermanyuk, Moulin, Pillet and Dauxois2017) cross-section, a vertical plate and a flat-topped hill (Brouzet et al. Reference Brouzet, Ermanyuk, Moulin, Pillet and Dauxois2017).
In nature, the main manifestation of monochromatic internal waves is the internal or baroclinic tide, generated in the ocean by the oscillation of the barotropic tide over bottom topography. The necessity of dealing with arbitrary topography has led to the development of a variety of theoretical approaches (Garrett & Kunze Reference Garrett and Kunze2007). Among them is the boundary integral method, and its numerical implementation the boundary element method. The forcing, typically a free-slip condition at a solid boundary, is replaced by a distribution of singularities on the boundary. The waves are expressed as the convolution of the distribution with the Green's function of the problem, and the boundary condition is reduced to an integral equation for the distribution. Once the distribution is known, the waves are obtained immediately in the whole fluid domain.
The boundary integral method is not specific to internal waves, or even to fluid mechanics. A review of its early history has been given by Cheng & Cheng (Reference Cheng and Cheng2005), and a recent review with extensive bibliography by Martin (Reference Martin2006, chapter 5). The numerical development of the method dates back to the 1960s, and the method is now used routinely for Stokes flow (Pozrikidis Reference Pozrikidis1992), potential flow (Pozrikidis Reference Pozrikidis2002) and acoustic waves (Crighton et al. Reference Crighton, Dowling, Ffowcs Williams, Heckl and Leppington1992, chapter 10), having reached the status of textbook topic (Hinch Reference Hinch2020, chapter 12). It is also commonly used in elasticity, heat transfer and electromagnetism; see, among many others, Brebbia, Telles & Wrobel (Reference Brebbia, Telles and Wrobel1984), Gaul, Kögl & Wagner (Reference Gaul, Kögl and Wagner2003) and Gibson (Reference Gibson2014).
For internal and inertial waves, the boundary integral method has been introduced first for the diffraction at a horizontal strip (Barcilon & Bleistein Reference Barcilon and Bleistein1969), a vertical strip (Robinson Reference Robinson1969), a wedge (Hurley Reference Hurley1970; Robinson Reference Robinson1970) and a shelf break (Hurley Reference Hurley1972). It was not pursued further in the Western literature at the time. In Russia, however, the derivation of Kirchhoff–Helmholtz integrals by Sobolev (Reference Sobolev1954) for inertial waves and Miropol'skii (Reference Miropol'skii1978) for internal waves, originally to solve the initial-value problem, led to investigation of the application of the same integrals to the boundary-value problem. This was done by Kapitonov (Reference Kapitonov1980) and Skazka (Reference Skazka1981) for inertial waves in three and two dimensions, respectively, Gabov & Shevtsov (Reference Gabov and Shevtsov1983, Reference Gabov and Shevtsov1984) for Boussinesq internal waves in three and two dimensions, respectively, and Gabov & Orazov (Reference Gabov and Orazov1986) for non-Boussinesq internal waves in one dimension, Pletner (Reference Pletner1991) and Sundukova (Reference Sundukova1991) in three dimensions and Pletner (Reference Pletner1992) and Allakhverdiev & Pletner (Reference Allakhverdiev and Pletner1993) in two dimensions.
A long series of papers followed, initiated by Gabov and continued by his collaborators after his death in 1989 (Sveshnikov, Shishmarev & Pletner Reference Sveshnikov, Shishmarev and Pletner1989). The series is presented in table 1, together with two separate papers by Korobkin (Reference Korobkin1990) and Davydova (Reference Davydova2006a), employing the same approach. The problem was considered for impulsive switch-on, and the steady state investigated in the large-time limit. Only the first papers in the series provided explicit solutions, while the later papers focused on the existence and unicity of a solution. The papers are largely derivative, exploring small variations of a limited number of configurations. Some results were published twice, in full form in regular journals and in summary form in Doklady. Some papers are identical (Kharik & Pletner Reference Kharik and Pletner1990a,Reference Kharik and Pletnerb), or only differ from each other by minute details (Krutitskii Reference Krutitskii1996a,Reference Krutitskiid, and also Reference Krutitskii1996e, Reference Krutitskii1997b), or adapt word-for-word an earlier study of inertia–gravity waves (Krutitskii Reference Krutitskii2000) to inertial waves (Krutitskii Reference Krutitskii2001) and internal waves (Krutitskii Reference Krutitskii2003b). In spite of these limitations, and because this body of work seems mostly unknown in the Western literature, it has felt useful to present it here.
Table 1. Publications based on the work of S.A. Gabov on the boundary integral method for internal and/or inertial waves, in two dimensions (2D) or three dimensions (3D). Unless stated otherwise (NB), the Boussinesq approximation is made throughout.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_tab1.png?pub-status=live)
Finally, three decades after its introduction for diffraction problems, the boundary integral method was brought back into the Western literature for generation problems, by Llewellyn Smith & Young (Reference Llewellyn Smith and Young2003) for the oscillations of a vertical barrier, considered analytically, and Pétrélis, Llewellyn Smith & Young (Reference Pétrélis, Llewellyn Smith and Young2006) for the oscillations of several topographies, considered numerically. Further analytical application was performed by Nycander (Reference Nycander2006) and Musgrave et al. (Reference Musgrave, Pinkel, MacKinnon, Mazloff and Young2016) to one or several vertical barriers, and further numerical application by Balmforth & Peacock (Reference Balmforth and Peacock2009), Echeverri & Peacock (Reference Echeverri and Peacock2010) and Echeverri et al. (Reference Echeverri, Yokossi, Balmforth and Peacock2011) to a variety of topographies. The scattering at a Gaussian topography was also considered numerically by Mathur, Carter & Peacock (Reference Mathur, Carter and Peacock2014). The mathematical foundations of the method, for both generation and scattering, were discussed by Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012).
In parallel, Sturova applied the boundary integral method to the oscillations of circular or elliptic horizontal cylinders, both analytically (Reference Sturova2001) and numerically (Reference Sturova2006, Reference Sturova2011), while Davydova & Chashechkin (Reference Davydova and Chashechkin2004) and Davydova (Reference Davydova2004, Reference Davydova2006b) discussed how the method can be combined with multiple scale analysis to calculate the waves and boundary layer generated by the oscillations of a vertical cylinder of arbitrary cross section in a slightly viscous fluid.
The aim of the present paper is twofold: first, to show how the method of coordinate stretching and analytic continuation can be combined with the boundary integral method, to solve the problem of internal wave generation by an oscillating body in an inviscid fluid; and second, to investigate the peculiarities of boundary integrals for monochromatic internal waves, compared with their usual properties for the Laplace and Helmholtz equations.
Boundary integrals have been considered by Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) already. We focus accordingly on two aspects complementary to their investigations: the equivalence (or lack thereof) between a direct formulation based on a Kirchhoff–Helmholtz integral involving both single and double layers, and indirect formulations based on single or double layers alone; and the discontinuities of the boundary integrals at the boundary.
The derivation of the waves is accompanied by that of a distribution of singularities equivalent to the body, with two consequences. First, once the spatial spectrum of the distribution is known, additional phenomena, like viscosity and unsteadiness, that affect the propagation of the waves and set their local structure, can be added into the analysis; this approach has been pioneered by Lighthill (Reference Lighthill1978, § 4.10) and Hurley & Keady (Reference Hurley and Keady1997), and recently developed by Voisin (Reference Voisin2020). Second, for a rigid body, the added mass of the body follows directly from the first moment of the distribution, providing immediate access to the radiated energy (called ‘conversion rate’ within the context of internal tides) and to the forces exerted on the body, without requiring the actual calculation of the waves; this aspect will be reported separately, and has been presented in summary form by Voisin (Reference Voisin2009).
The problem of internal wave generation by an oscillating body is stated in § 2 and its direct formulation as a boundary integral in § 3. It is solved in §§ 4 and 5 for an elliptic cylinder and a spheroid, as prototypal two- and three-dimensional bodies, respectively. Arbitrary oscillations are considered first, and the results applied to the two simplest types of oscillations: radial pulsations and rigid vibrations. The indirect formulation of the problem as a single-layer integral is presented in § 6 and shown to lead to a simple representation of the body. The representations of a vibrating sphere in Voisin, Ermanyuk & Flór (Reference Voisin, Ermanyuk and Flór2011) and a vibrating elliptic cylinder in Voisin (Reference Voisin2020) are recovered as particular cases. Sections 7 and 8 discuss the behaviour of the boundary integrals at the boundary and in the far field, respectively. They are followed by a conclusion in § 9.
2. Statement of the problem
We consider linear internal waves in an inviscid uniformly stratified Boussinesq fluid of buoyancy frequency $N$, and use the generalized potential introduced by Sobolev (Reference Sobolev1954) for inertial waves and Gorodtsov & Teodorovich (Reference Gorodtsov and Teodorovich1980), Hart (Reference Hart1981), Gray, Hart & Farrell (Reference Gray, Hart and Farrell1983), Gabov & Mamedov (Reference Gabov and Mamedov1983), Gabov, Malysheva & Sveshnikov (Reference Gabov, Malysheva and Sveshnikov1983), Gabov et al. (Reference Gabov, Malysheva, Sveshnikov and Shatov1984), Voisin (Reference Voisin1991, Reference Voisin2003) and Gorodtsov (Reference Gorodtsov2013) for internal waves, with specifics listed in table 2. Gauge invariance, ensuring the completeness of the representation, has been discussed by Sobolev (Reference Sobolev1954), Gray et al. (Reference Gray, Hart and Farrell1983) and Gabov et al. (Reference Gabov, Malysheva and Sveshnikov1983). A different but related representation based on the toroidal/poloidal decomposition has been introduced by Kistovich & Chashechkin (Reference Kistovich and Chashechkin2001). Denoting the potential, called ‘internal’ by Voisin (Reference Voisin1991), as
$\psi$, the wave equation takes the form
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn1.png?pub-status=live)
with $z$ the vertical coordinate,
$\boldsymbol {x} = (x,y,z)$ the position,
$\boldsymbol {\nabla } = (\partial /\partial x,\partial /\partial y,\partial /\partial z)$ the del operator and
$\boldsymbol {\nabla }_{{h}} = (\partial /\partial x,\partial /\partial y,0)$ its horizontal projection. The forcing is represented as a source of mass releasing the volume
$q = \boldsymbol {\nabla }\boldsymbol {\cdot }\boldsymbol {u}$ of fluid per unit volume per unit time. The disturbances
$\boldsymbol {u}$ in velocity,
$p$ in pressure and
$\rho$ in density are expressed in terms of
$\psi$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn2.png?pub-status=live)
with $\rho _0$ the density at rest and
$g$ the acceleration due to gravity.
Table 2. Introduction of generalized potentials for waves in rotating and/or stratified fluids with or without viscosity, non-Boussinesq effects, compressibility and heat conduction.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_tab2.png?pub-status=live)
We are interested in monochromatic waves of frequency $\omega$, varying in time through the factor
$\exp (-\mathrm {i}\omega t)$ which is suppressed in the following. The wave equation becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn3.png?pub-status=live)
and the fluid dynamical quantities become
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn4.png?pub-status=live)
An oscillating body of surface $S$ and outward normal
$\boldsymbol {n}$ generates waves by imposing a normal velocity
$U_n$ at
$S$, through the free-slip boundary condition
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn5.png?pub-status=live)
The boundary integral method replaces this condition by a source term $q$ in the wave equation.
3. Direct approach
3.1. Kirchhoff–Helmholtz integral
For this, we adapt the approach used for the Laplace and Helmholtz equations; see, for example, Lighthill (Reference Lighthill1986, § 8.1), Jackson (Reference Jackson1999, § 1.8), Martin (Reference Martin2006, § 5.6) and Pierce (Reference Pierce2019, § 4.6). For two arbitrary functions $f$ and
$g$ we write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn6.png?pub-status=live)
which upon integration inside a volume $V$ delimited by the surface
$S$ of outward normal
$\boldsymbol {n}$, and application of the divergence theorem, gives
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn7.png?pub-status=live)
where $\partial /\partial n = \boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {\nabla }$ and
$\partial /\partial n_{{h}} = \boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {\nabla }_{{h}}$.
We apply this Green's theorem to the volume $V_+$ exterior to the oscillating body, delimited internally by the surface
$S$ of the body (with outward normal
$-\boldsymbol {n}$) and externally by a surface
$S_\infty$ at infinity, as shown in figure 1. The volume interior to the body is denoted as
$V_-$. We choose a fixed observation point
$\boldsymbol {x}$ in either
$V_+$ or
$V_-$, and denote the variable integration point in
$V_+$ as
$\boldsymbol {x}'$. We take for
$f$ the internal potential
$\psi (\boldsymbol {x}')$ at
$\boldsymbol {x}'$, and for
$g$ the propagator
$G(\boldsymbol {x}-\boldsymbol {x}')$ from
$\boldsymbol {x}'$ to
$\boldsymbol {x}$. Here,
$G$ is the Green's function of the wave equation, corresponding to unit point forcing and satisfying
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn8.png?pub-status=live)
with $\delta$ the Dirac delta function, together with the causality condition that
$G$ be analytic in the upper half of the complex
$\omega$-plane. We assume that both
$G$ and
$\psi$ decay fast enough at infinity for the contribution of
$S_\infty$ to vanish, and will come back to this assumption later in § 8. We obtain the Kirchhoff–Helmholtz integral
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn9.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn10.png?pub-status=live)
where $\partial /\partial n' = \boldsymbol {n}'\boldsymbol {\cdot }\boldsymbol {\nabla }'$ and
$\partial /\partial n'_{{h}} = \boldsymbol {n}'\boldsymbol {\cdot }\boldsymbol {\nabla }'_{{h}}$, with
$\boldsymbol {n}' = \boldsymbol {n}(\boldsymbol {x}')$.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_fig1.png?pub-status=live)
Figure 1. Surfaces for the derivation of the Kirchhoff–Helmholtz integral.
Letting now $\boldsymbol {x}$ approach
$S$ from within
$V_+$, we obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn11.png?pub-status=live)
an integral relation between the surface values $\varPsi$ of the internal potential (hence the pressure) and
$U_n$ of the normal velocity, to be satisfied for all
$\boldsymbol {x}$ on the outer side
$S_+$ of
$S$. Given the pressure this provides an integral equation of the first kind for the velocity, and given the velocity an equation of the second kind for the pressure. Once it is solved, the waves are expressed through the fluid, for
$\boldsymbol {x} \in V_+$, as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn12.png?pub-status=live)
This expression combines two terms: a single-layer potential
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn13.png?pub-status=live)
generated by the surface distribution of monopoles
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn14.png?pub-status=live)
with density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn15.png?pub-status=live)
where $\delta _S$ is the Dirac delta function of support
$S$, such that, if
$S(\boldsymbol {x}) = 0$ is the equation of the surface,
$\delta _S(\boldsymbol {x}) = |\boldsymbol {\nabla } S|\delta [S(\boldsymbol {x})]$; and a double-layer potential
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn16.png?pub-status=live)
generated by the surface distribution of dipoles
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn17.png?pub-status=live)
with density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn18.png?pub-status=live)
See, for example, Jackson (Reference Jackson1999, § 1.6) for the introduction of such potentials in electrostatics, and Martin (Reference Martin2006, § 5.3) for non-dispersive waves.
The forcing is thus represented by the equivalent source
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn19.png?pub-status=live)
Such representation allows the extension of (3.6) to a viscous fluid, once the spectrum
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn20.png?pub-status=live)
is known, namely
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn21.png?pub-status=live)
with $\boldsymbol {k}$ the wave vector and
$\boldsymbol {k}_{{h}}$ its horizontal projection. The extension follows the lines laid by Lighthill (Reference Lighthill1978, § 4.10), Hurley & Keady (Reference Hurley and Keady1997) and Voisin (Reference Voisin2020); it will be discussed further in § 6.2.
3.2. Green's function
To determine the Green's function we use the technique introduced by Bryan (Reference Bryan1889) for inertial waves and Hurley (Reference Hurley1972) for internal waves. In the frequency range $\omega > N$ of evanescent waves, (3.3), rewritten as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn22.png?pub-status=live)
is elliptic. Stretching the coordinates according to
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn23.png?pub-status=live)
transforms it into a Poisson equation, of known Green's function (Bleistein Reference Bleistein1984, § 6.2). Causality allows the continuation of the solution to the frequency range $0 < \omega < N$ of propagating waves. Continuation is implemented via Reference LighthillLighthill's (Reference Lighthill1978, § 4.9) radiation condition, namely by adding to the frequency a small positive imaginary part
$\epsilon$ which is later allowed to tend to
$0$; in other words, by performing the replacement
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn24.png?pub-status=live)
The Green's function is obtained in three dimensions as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn25.png?pub-status=live)
that is in unstretched coordinates
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn26.png?pub-status=live)
with associated velocity
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn27.png?pub-status=live)
Similarly, in two dimensions, when the waves are independent of the horizontal coordinate $y$, the Green's function is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn28.png?pub-status=live)
Knowing the Green's function we can now solve (3.5) for prototypal oscillating bodies, namely an elliptic cylinder in two dimensions and a spheroid in three dimensions, in the next two sections.
4. Elliptic cylinder
An elliptic cylinder of horizontal semi-axis $a$, vertical semi-axis
$b$ and equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn29.png?pub-status=live)
is transformed by the stretching (3.17) into another elliptic cylinder, of semi-axes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn30.png?pub-status=live)
respectively. For $(a/b)^2+(N/\omega )^2 > 1$, corresponding to an original ellipse which either has its major axis horizontal and operates at any frequency
$\omega > N$, or has its major axis vertical and operates in the frequency range
$N < \omega < N/[1-(a/b)^2]^{1/2}$, the stretched ellipse has its major axis horizontal. We consider this situation first, and will deal with the other situations by analytic continuation.
4.1. Solution for
$\omega > N$ and
$(a/b)^2+(N/\omega )^2 > 1$
For $(a/b)^2+(N/\omega )^2 > 1$, the stretched ellipse has its foci along the horizontal axis, at
$(x_\star = \pm c, z_\star = 0)$, with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn31.png?pub-status=live)
We introduce elliptic coordinates $(\xi ,\eta$) defined by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn32.png?pub-status=live)
with $0 \leqslant \xi < \infty$ and
$0 \leqslant \eta < 2{\rm \pi}$. This definition is inverted in terms of the distances to the foci,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn33.png?pub-status=live)
as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn34.png?pub-status=live)
see Happel & Brenner (Reference Happel and Brenner1983, Appendix A) or Landau & Lifshitz (Reference Landau and Lifshitz1984, § 4). The spatial derivatives become
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn35.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn36.png?pub-status=live)
where $\sinh ^2\xi +\sin ^2\eta = \sinh (\xi +\mathrm {i}\eta )\sinh (\xi -\mathrm {i}\eta )$.
At the ellipse, $\xi = \xi _0$ with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn37.png?pub-status=live)
so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn38.png?pub-status=live)
There, $\eta$ reduces to the eccentric angle for the original ellipse, such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn39.png?pub-status=live)
see Sommerville (Reference Sommerville1933, § IV.10) or Milne-Thomson (Reference Milne-Thomson1968, § 6.32). The arc length element follows as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn40.png?pub-status=live)
and the outward normal as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn41.png?pub-status=live)
where $\boldsymbol {e}_x$ and
$\boldsymbol {e}_z$ are unit vectors along the
$x$- and
$z$-axes, respectively. The boundary integral equation (3.5) becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn42.png?pub-status=live)
where $\xi = \xi _0+0 > \xi ' = \xi _0$.
To solve it, we adapt the procedure used by Gorodtsov & Teodorovich (Reference Gorodtsov and Teodorovich1982) for the potential flow around a circular cylinder. The normal velocity $U_n(\eta )$ is expanded in circular functions as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn43.png?pub-status=live)
with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn44.png?pub-status=live)
and the surface potential $\varPsi (\eta )$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn45.png?pub-status=live)
with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn46.png?pub-status=live)
Here, $\epsilon _m = 1$ for
$m = 0$ and
$2$ for
$m \geqslant 1$ is the Neumann factor. The Green's function has been shown by Morse & Feshbach (Reference Morse and Feshbach1953, p. 1202) to admit of the expansion
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn47.png?pub-status=live)
where $\xi _< = \min (\xi ,\xi ')$ and
$\xi _> = \max (\xi ,\xi ')$. Use of these expansions turns (4.13) into a diagonal (infinite) linear system, of solution
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn48.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn49.png?pub-status=live)
Evaluation of the convolution integral (3.6), followed by differentiation according to (2.4), yields
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn50.png?pub-status=live)
for the pressure, and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn51.png?pub-status=live)
for the velocity.
4.2. Solution for
$\omega > N$ and
$(a/b)^2+(N/\omega )^2 < 1$
These results are continued analytically to $(a/b)^2+(N/\omega )^2 < 1$ by applying (3.18). The stretched ellipse has its foci along the vertical axis, at
$(x_\star = 0,z_\star = \pm c')$, with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn52.png?pub-status=live)
This prompts the introduction of elliptic coordinates $(\xi ',\eta )$, defined by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn53.png?pub-status=live)
where $0 \leqslant \xi ' < \infty$ and
$0 \leqslant \eta < 2{\rm \pi}$. Writing
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn54.png?pub-status=live)
the solution of the problem is seen to follow immediately from that in § 4.1, by replacing $\xi$ by
$\xi '+\mathrm {i}{\rm \pi} /2$ and
$c$ by
$-\mathrm {i}c'$ throughout.
4.3. Solution for
$\omega < N$
In the frequency range $0 < \omega < N$, the waves propagate at the angle
$\theta _0 = \arccos (\omega /N)$ to the vertical. The semi-focal distance of the stretched ellipse becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn55.png?pub-status=live)
a quantity interpreted by Hurley (Reference Hurley1997) as the half-width of the wave beams delimited by the critical wave rays tangential to the ellipse on either side. The semi-axes become
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn56.png?pub-status=live)
where
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn57.png?pub-status=live)
is the eccentric angle of the critical points at which the critical rays are tangential to the ellipse. Both definitions are illustrated in figure 2.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_fig2.png?pub-status=live)
Figure 2. Geometry for the beams of waves propagating (a) upward to the right and downward to the left, and (b) downward to the right and upward to the left.
The stretched elliptic coordinates $\xi$ and
$\eta$ become complex, with
$\xi$ taking the value
$\xi _0 = \mathrm {i}\eta _0$ at the ellipse. The greatest difficulty here is the expression of
$\xi$ and
$\eta$ in more comprehensible coordinates. For a sphere or spheroid, Sarma & Krishna (Reference Sarma and Krishna1972) expressed the associated stretched spheroidal coordinates in terms of algebraic functions of the Cartesian coordinates, but did not consider the determination of these multivalued functions explicitly. Hendershott (Reference Hendershott1969) and Voisin (Reference Voisin1991) set the determination in the far field. Appleby & Crighton proceeded differently, for both the circular cylinder (Reference Appleby and Crighton1986) and the sphere (Reference Appleby and Crighton1987), decomposing the wave field into zones delimited by the critical rays, and investigating the properties of the stretched coordinates in each zone, concluding for the sphere that ‘this illustrates the difficulty of working in more comprehensible coordinates’. Davis (Reference Davis2012) combined the two approaches together for the sphere, expressing, in each zone, the stretched coordinates in terms of real algebraic functions of the Cartesian coordinates.
In actuality, Lighthill's radiation condition (3.18) allows both the stretched coordinates to be expressed in terms of the characteristic coordinates
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn58.png?pub-status=live)
shown in figure 2, and the associated multivalued functions to be determined unequivocally. The inversion formulae (4.5)–(4.6) are not convenient for this purpose. We proceed instead from (3.17), (4.3) and (4.4), writing
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn59.png?pub-status=live)
Expanding these for $\omega = N\cos \theta _0+\mathrm {i}\epsilon$, with
$0 < \epsilon /N \ll 1$, we get
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn60.png?pub-status=live)
so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn61.png?pub-status=live)
In addition to the quantities associated with the wave beams, namely their angle $\theta _0$ to the vertical, their half-width
$c$ and the coordinates
$(x_\pm ,z_\pm )$ perpendicular to and along them, respectively, this expression involves also quantities associated with the critical segments joining, for each beam, the opposite critical points on either side of the ellipse, namely their angle
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn62.png?pub-status=live)
to the horizontal, their half-length
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn63.png?pub-status=live)
and the coordinates
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn64.png?pub-status=live)
along and perpendicular to them, respectively. These quantities, introduced by Voisin (Reference Voisin2020), are represented in figure 2.
We may then write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn65.png?pub-status=live)
on the understanding that the determination of the square roots is set by the replacement
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn66.png?pub-status=live)
so that, in particular,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn67.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn68.png?pub-status=live)
This determination, shown in figure 3, coincides with those in figures 3–4 of Hurley (Reference Hurley1972) for a circular cylinder and figure 3 of Hurley (Reference Hurley1997) for an elliptic cylinder.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_fig3.png?pub-status=live)
Figure 3. (a) Determinations of $(x_+^2-c^2)^{1/2}$ and (b) determinations of
$(x_-^2-c^2)^{1/2}$.
It thus follows that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn69.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn70.png?pub-status=live)
leading for the pressure to the expansion
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn71.png?pub-status=live)
and for the velocity to
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn72.png?pub-status=live)
with $\boldsymbol {e}_{z_\pm } = \pm \boldsymbol {e}_x\sin \theta _0+\boldsymbol {e}_z\cos \theta _0$ a unit vector along the
$z_\pm$-axis. Both expansions are of the form anticipated by Barcilon & Bleistein (Reference Barcilon and Bleistein1969) and Hurley (Reference Hurley1972) for a circular cylinder.
At the ellipse (of contour $C$ with outer side
$C_+$), the pressure becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn73.png?pub-status=live)
and the velocity becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn74.png?pub-status=live)
Both the pressure and the normal velocity are regular, the latter being equal to its prescribed value $U_n$. The tangential velocity is singular at the critical points
$\eta = \eta _0$,
${\rm \pi} -\eta _0$,
${\rm \pi} +\eta _0$ and
$2{\rm \pi} -\eta _0$; there, even at vanishingly small viscosity, the actual no-slip boundary condition will come into play, giving rise to the boundary-layer eruption predicted by Kerswell (Reference Kerswell1995) and Le Dizès & Le Bars (Reference Le Dizès and Le Bars2017).
4.4. Particular cases
We consider now, as is usually done in acoustics, the simplest types of oscillations, associated with the lowest values of $m$; see Lighthill (Reference Lighthill1978, § 1.11) or Pierce (Reference Pierce2019, §§ 4.1–2).
Monopolar oscillations, for which $m = 0$, correspond to radial pulsations at the velocity
$\boldsymbol {U} = U\boldsymbol {x}/c$, such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn75.png?pub-status=live)
The associated pressure is uniform at the cylinder, with value
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn76.png?pub-status=live)
This mode of oscillation, which in the presence of viscosity gives the lowest rate of decrease of the wave amplitude with distance away from the cylinder, is the least studied in the laboratory, having only been considered by Makarov, Neklyudov & Chashechkin (Reference Makarov, Neklyudov and Chashechkin1990) and Machicoane et al. (Reference Machicoane, Cortet, Voisin and Moisy2015) for a circular cylinder.
Dipolar oscillations, for which $m = 1$, correspond to back-and-forth vibrations of a rigid cylinder at the velocity
$\boldsymbol {U} = U\boldsymbol {e}_x+W\boldsymbol {e}_z$, such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn77.png?pub-status=live)
The pressure varies linearly with the Cartesian coordinates at the cylinder, in the form
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn78.png?pub-status=live)
This configuration is the most studied in the laboratory, from the visualizations of Mowbray & Rarity (Reference Mowbray and Rarity1967) up to the quantitative measurements of Sutherland et al. (Reference Sutherland, Dalziel, Hughes and Linden1999) and Zhang, King & Swinney (Reference Zhang, King and Swinney2007) for a circular cylinder and Sutherland & Linden (Reference Sutherland and Linden2002) for elliptic cylinders of aspect ratios $a/b = 1$,
$2$ and
$3$.
The expressions of the pressure and velocity at propagating frequencies $0 < \omega < N$ are given in table 3. For the vibrating cylinder they involve the notation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn79.png?pub-status=live)
introduced by Hurley (Reference Hurley1997). They are accompanied by the limits $a \to 0$, corresponding to the horizontal vibrations of a vertical knife edge, considered theoretically by Llewellyn Smith & Young (Reference Llewellyn Smith and Young2003) and experimentally by Peacock, Echeverri & Balmforth (Reference Peacock, Echeverri and Balmforth2008), and
$b \to 0$, corresponding to the vertical vibrations of a horizontal knife edge.
Table 3. Pressure and velocity in the propagating frequency range $0 < \omega < N$. The determination of the square roots is set according to (4.36), becoming
$x_\pm \to x_\pm +\mathrm {i}\operatorname {sign} x$ and
$x_\pm \to x_\pm \pm \mathrm {i}\operatorname {sign} z$ for vertical and horizontal knife edges, respectively.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_tab3.png?pub-status=live)
5. Spheroid
We move on to a spheroid of horizontal semi-axis $a$, vertical semi-axis
$b$ and equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn80.png?pub-status=live)
The approach remains the same as for the elliptic cylinder, but the exposition is made more intricate by the switch from hyperbolic functions to Legendre functions, and from circular functions to spherical harmonics. The relevant properties of these functions are recalled in Appendices A and B; they will be used silently in the following.
5.1. Solution for
$\omega > N$ and
$(a/b)^2+(N/\omega )^2 > 1$
For $(a/b)^2+(N/\omega )^2 > 1$, the coordinate stretching (3.17) transforms the original spheroid into an oblate one, of semi-axes
$a_\star$ and
$b_\star$ given by (4.2), and focal circle
$(r_\star = c,z_\star = 0)$. Here,
$r_{{h}} = (x^2+y^2)^{1/2}$ is the horizontal radial distance, stretched as
$r_\star = (x_\star ^2+y_\star ^2)^{1/2}$, and
$c$ is given by (4.3). We introduce oblate spheroidal coordinates
$(\xi ,\eta ,\phi )$ defined by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn81.png?pub-status=live)
where $0 \leqslant \xi < \infty$ and
$0 \leqslant \eta \leqslant {\rm \pi}$, with
$0 \leqslant \phi < 2{\rm \pi}$ the usual azimuthal angle, such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn82.png?pub-status=live)
In a given azimuthal plane, this definition is inverted in terms of the distances to the two foci in this plane,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn83.png?pub-status=live)
as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn84.png?pub-status=live)
see Happel & Brenner (Reference Happel and Brenner1983, Appendix A) or Landau & Lifshitz (Reference Landau and Lifshitz1984, § 4). Differentiation is expressed as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn85.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn86.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn87.png?pub-status=live)
where $\sinh ^2\xi +\cos ^2\eta = \cosh (\xi +\mathrm {i}\eta )\cosh (\xi -\mathrm {i}\eta )$.
At the spheroid, $\xi = \xi _0$ with
$\xi _0$ given by (4.9), and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn88.png?pub-status=live)
implying that ${\rm \pi} /2-\eta$ is Legendre's (Reference Legendre1806) reduced latitude and Cayley's (Reference Cayley1870) parametric latitude. The surface area element follows as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn89.png?pub-status=live)
with $\mathrm {d}\varOmega = \sin \eta \,\mathrm {d}\eta \,\mathrm {d}\phi$ the solid angle element, and the outward normal as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn90.png?pub-status=live)
where $\boldsymbol {e}_{r_{{h}}} = \boldsymbol {e}_x\cos \phi +\boldsymbol {e}_y\sin \phi$ and
$\boldsymbol {e}_\phi = -\boldsymbol {e}_x\sin \phi +\boldsymbol {e}_y\cos \phi$ are unit vectors along the radial and azimuthal horizontal directions, respectively, with
$\boldsymbol {e}_y$ a unit vector along the
$y$-axis. The boundary integral equation (3.5) becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn91.png?pub-status=live)
where $\xi = \xi _0+0 > \xi ' = \xi _0$.
To solve it, we adapt the procedure used by Gorodtsov & Teodorovich (Reference Gorodtsov and Teodorovich1982) for the potential flow around a sphere. The normal velocity $U_n(\eta ,\phi )$ is expanded in spherical harmonics as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn92.png?pub-status=live)
with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn93.png?pub-status=live)
where an overbar denotes a complex conjugate, and the surface potential $\varPsi (\eta ,\phi )$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn94.png?pub-status=live)
with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn95.png?pub-status=live)
The expansion of the Green's function has been given by Morse & Feshbach (Reference Morse and Feshbach1953, p. 1296) and Hobson (Reference Hobson1931, § 251), in both cases with typos: an extraneous factor $2$ for the former, and a missing factor
$\mathrm {i}$ for the latter. Accordingly, the expansion has been rederived by an adaptation of the procedure of Jackson (Reference Jackson1999, § 3.9), to get
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn96.png?pub-status=live)
where $\xi _< = \min (\xi ,\xi ')$ and
$\xi _> = \max (\xi ,\xi ')$. The solution of (5.10) is then straightforward, in the form
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn97.png?pub-status=live)
Convolution according to (3.6), followed by differentiation according to (2.4), yields
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn98.png?pub-status=live)
for the pressure, and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn99.png?pub-status=live)
for the velocity.
5.2. Solution for
$\omega > N$ and
$(a/b)^2+(N/\omega )^2 < 1$
For $(a/b)^2+(N/\omega )^2 < 1$, we proceed as in § 4.2. The stretched spheroid becomes prolate, having foci at
$(r_\star = 0,z_\star = \pm c')$, with
$c'$ given by (4.22). This leads to the introduction of prolate spheroidal coordinates
$(\xi ',\eta ,\phi )$ defined by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn100.png?pub-status=live)
with $0 \leqslant \xi ' < \infty$ and
$0 \leqslant \eta \leqslant {\rm \pi}$, so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn101.png?pub-status=live)
Accordingly, the waves follow from those in § 5.1, by replacing $\xi$ by
$\xi '+\mathrm {i}{\rm \pi} /2$ and
$c$ by
$-\mathrm {i}c'$ throughout.
5.3. Solution for
$\omega < N$
For $0 < \omega < N$, the waves propagate in beams inclined at the angle
$\theta _0 = \arccos (\omega /N)$ to the vertical. The beams have half-width
$c$ given by (4.25), and are delimited by the critical rays tangential to the spheroid above and below, forming two double cones grazing the spheroid at critical circles of reduced latitude
$\eta _0$ given by (4.27). The spheroid becomes the surface
$\xi = \mathrm {i}\eta _0$, and the waves are deduced by applying the transformation (3.18) to (5.17)–(5.18), writing
$\xi _0 \to \mathrm {i}\eta _0+0$ so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn102.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn103.png?pub-status=live)
For the complex coordinates $\xi$ and
$\eta$ we proceed as in § 4.3, writing
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn104.png?pub-status=live)
and thence
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn105.png?pub-status=live)
where
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn106.png?pub-status=live)
are the characteristic coordinates associated with the wave beams, and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn107.png?pub-status=live)
are the coordinates associated with the truncated double cone joining the two critical circles at the surface of the spheroid. Here, $\chi _0$, given by (4.32), is the angle of the generatrices of the double cone to the horizontal (in other words, the critical latitude), and
$d$, given by (4.33), their half-length.
We then write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn108.png?pub-status=live)
where the determination of the square roots is set by the replacement (4.36), yielding the combinations
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn109.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn110.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn111.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn112.png?pub-status=live)
consistent with the determinations in § 4 of Davis (Reference Davis2012). The pressure is obtained as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn113.png?pub-status=live)
and the velocity as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn114.png?pub-status=live)
where $\boldsymbol {e}_{z_\pm } = \pm \boldsymbol {e}_{r_{{h}}}\sin \theta _0+\boldsymbol {e}_z\cos \theta _0$ is a unit vector along the
$z_\pm$-axis. The three- dimensional geometry breaks the relative simplicity of two-dimensional waves: the waves can no longer be separated into two components, one depending exclusively on
$x_+$ and the other on
$x_-$; instead, they are expressed in terms of the combinations (5.28)–(5.31), each of which involves both
$x_+$ and
$x_-$.
5.4. Particular cases
The simplest mode of oscillation of the spheroid is radial pulsations at the velocity $\boldsymbol {U} = U\boldsymbol {x}/c$. Such monopolar forcing is of degree
$l = 0$, with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn115.png?pub-status=live)
The expressions of the pressure and velocity at propagating frequencies $0 < \omega < N$ are given in table 4. They are consistent with Hendershott (Reference Hendershott1969), Appleby & Crighton (Reference Appleby and Crighton1987), Voisin (Reference Voisin1991), Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) and Davis (Reference Davis2012) for a sphere. The pressure is uniform at the spheroid,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn116.png?pub-status=live)
as is expected for pulsations.
Table 4. Pressure and velocity in the propagating frequency range $0 < \omega < N$. The determination of the square roots is set according to (4.36), becoming
$x_\pm \to x_\pm \pm \mathrm {i}\operatorname {sign} z$ for the horizontal disc.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_tab4.png?pub-status=live)
When the spheroid is rigid and vibrates back-and-forth with velocity $\boldsymbol {U} = U\boldsymbol {e}_x+V\boldsymbol {e}_y+W\boldsymbol {e}_z$, the forcing is dipolar and of degree
$l = 1$, with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn117.png?pub-status=live)
Introducing the ratio $\varUpsilon = \tanh \xi _0$, that is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn118.png?pub-status=live)
and the combination $D(\varUpsilon ) = \cosh ^2\xi _0(1-\sinh \xi _0\operatorname {arccot}\sinh \xi _0)$, that is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn119.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn120.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn121.png?pub-status=live)
the pressure and velocity become as given in table 4 for $0 < \omega < N$. They are consistent with Sarma & Krishna (Reference Sarma and Krishna1972) and Lai & Lee (Reference Lai and Lee1981) for a spheroid, and Appleby & Crighton (Reference Appleby and Crighton1987), Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) and Davis (Reference Davis2012) for a sphere. Their relation to the experiments of Flynn, Onu & Sutherland (Reference Flynn, Onu and Sutherland2003), King, Zhang & Swinney (Reference King, Zhang and Swinney2009), Voisin et al. (Reference Voisin, Ermanyuk and Flór2011) and Ghaemsaidi & Peacock (Reference Ghaemsaidi and Peacock2013) for a sphere will be discussed in § 6.2. At the spheroid, the pressure exhibits linear variations with the Cartesian coordinates,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn122.png?pub-status=live)
The limit $a \to 0$ of the vibrating spheroid corresponds to a vertical needle, which generates no waves, and the limit
$b \to 0$ to the heaving oscillations of a horizontal circular disc, considered theoretically by Sarma & Krishna (Reference Sarma and Krishna1972), Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2011) and Davis (Reference Davis2012) in an inviscid fluid and Davis & Llewellyn Smith (Reference Davis and Llewellyn Smith2010) in a viscous fluid, and experimentally by Bardakov, Vasil'ev & Chashechkin (Reference Bardakov, Vasil'ev and Chashechkin2007). For the disc, we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn123.png?pub-status=live)
yielding the pressure and velocity in table 4.
6. Indirect approach
6.1. Layer potentials
Investigations so far have been based on the Kirchhoff–Helmholtz integral (3.4), combining the single-layer potential (3.7) and the double-layer potential (3.10). All the quantities involved are physically meaningful, either internal potential (hence pressure) inside the fluid in $V_+$, or internal potential and normal velocity on the side
$S_+$ of the surface
$S$ of the body in contact with the fluid. By contrast, the layer potentials are purely mathematical entities, defined not only in the fluid but also inside the oscillating body in
$V_-$. For each potential, the associated source density involves values of the potential and its derivative on both sides of
$S$, inner for
$S_-$ and outer for
$S_+$.
Specifically, the single-layer potential satisfies the differential equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn124.png?pub-status=live)
implying (Schwartz Reference Schwartz1966, § II.2.3) continuity of the potential and discontinuity of its modified normal derivative at $S$, namely
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn125.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn126.png?pub-status=live)
As a consequence, the monopole density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn127.png?pub-status=live)
is equal to the discontinuity of the normal velocity at $S$, going from
$S_-$ to
$S_+$. Similarly, the double-layer potential satisfies
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn128.png?pub-status=live)
implying its discontinuity and the continuity of its modified normal derivative at $S$, with
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn129.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn130.png?pub-status=live)
As a consequence, the dipole density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn131.png?pub-status=live)
is opposite to the discontinuity of the internal potential at $S$, going from
$S_-$ to
$S_+$.
By imagining the oscillating body to be filled with fluid, it is possible to give physical meaning to each layer potential and to consider representing the motion of the whole fluid by this potential alone, without the other. For this, we adapt the procedure described for the Laplace and Helmholtz equations by Lamb (Reference Lamb1932, §§ 57, 58 and 290) and Copley (Reference Copley1968). The volume $V_-$ of the body is replaced by fictitious stratified fluid, of the same buoyancy frequency
$N$ as the real fluid outside. Applying Green's theorem (3.2) inside
$V_-$, we obtain for the internal potential
$\psi _{{f}}$ the Kirchhoff–Helmholtz integral
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn132.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn133.png?pub-status=live)
which combined with (3.4) gives
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn134.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn135.png?pub-status=live)
The imposition of the boundary condition $\psi _{{f}} = \psi$ at
$S$ yields a single-layer field of monopole density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn136.png?pub-status=live)
and similarly the condition $\boldsymbol {u}_{{f}}\boldsymbol {\cdot }\boldsymbol {n} = \boldsymbol {u}\boldsymbol {\cdot }\boldsymbol {n}$ yields a double-layer field of dipole density
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn137.png?pub-status=live)
In both cases the wave field is described by $\psi$ outside
$S$ and
$\psi _{{f}}$ inside it.
A calculation of the waves based on the Kirchhoff–Helmholtz integral (3.4) is called direct, and a calculation based on the single-layer potential (3.7) or double-layer potential (3.10) is called indirect. The legitimacy of direct approaches follows from the derivation of the Kirchhoff–Helmholtz integral in § 3.1, but the legitimacy of indirect approaches remains to be seen. Paraphrasing Lamb (Reference Lamb1932, § 290), it would be wrong to assume that, as in the case of the ordinary potential, a function $\psi _{{f}}$ necessarily exists which satisfies the wave equation (2.3) throughout the finite region
$V_-$, and also fulfils the condition that
$\psi _{{f}}$ or its modified normal derivative shall assume arbitrarily prescribed values over the boundary
$S$.
6.2. Equivalent source
We consider the single layer (3.7), associated with the equivalent source (3.8). The boundary condition (2.5) yields for this source the integro-differential equation
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn138.png?pub-status=live)
to be satisfied for all $\boldsymbol {x}\in S_+$.
In two dimensions, for the elliptic cylinder of § 4, the equation becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn139.png?pub-status=live)
where $\xi = \xi _0+0 > \xi ' = \xi _0$. Expanding the source density
$\sigma (\eta )$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn140.png?pub-status=live)
and the normal velocity $U_n(\eta )$ as (4.14), the solution is obtained immediately as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn141.png?pub-status=live)
The waves (4.20)–(4.21) then follow from (3.7) and (2.4).
The cylinder is thus equivalent to the source
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn142.png?pub-status=live)
of spectrum
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn143.png?pub-status=live)
where $\operatorname {J}_m$ denotes a cylindrical Bessel function,
$\boldsymbol {k} = (k_x,k_z)$ the wave vector and
$\eta _k$ the angle such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn144.png?pub-status=live)
The knowledge of this spectrum provides an alternative expression of the waves, into which the effect of viscous attenuation can be added following Lighthill (Reference Lighthill1978, § 4.10) at large distance from the source, and Voisin (Reference Voisin2020) at arbitrary distance from it.
We use the latter. In an inviscid fluid, the velocity outside a source of mass $q(x,z)$ of support the elliptic domain
$x^2/a^2+z^2/b^2 < 1$ is given in the propagating frequency range
$0 < \omega < N$ by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn145.png?pub-status=live)
where
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn146.png?pub-status=live)
and $q_\pm (k_{x_\pm },k_{z_\pm }) = q(k_x,k_z)$. For the cylinder, the spectrum (6.16) becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn147.png?pub-status=live)
so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn148.png?pub-status=live)
Using the transform
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn149.png?pub-status=live)
taken from table 5 of Voisin (Reference Voisin2003), the velocity (4.41) is recovered. Viscosity, at large Stokes number $\omega c^2/\nu \gg 1$ with
$\nu$ the kinematic viscosity, adds an exponential attenuation factor inside the integrand of (6.21), transforming it into
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn150.png?pub-status=live)
with $\beta = \nu /(2N\sin \theta _0)$.
In three dimensions, for the spheroid of § 5, the integro-differential equation (6.11) becomes
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn151.png?pub-status=live)
where $\xi = \xi _0+0 > \xi ' = \xi _0$. Expanding
$\sigma (\eta ,\phi )$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn152.png?pub-status=live)
and $U_n(\eta ,\phi )$ as (5.11), we obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn153.png?pub-status=live)
which leads to the waves (5.17)–(5.18). The spheroid is thus represented by the source
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn154.png?pub-status=live)
Its spectrum is obtained using (A14) as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn155.png?pub-status=live)
where $\operatorname {j}_l$ denotes a spherical Bessel function,
$\boldsymbol {k} = (k_x,k_y,k_z)$ the wave vector, of horizontal modulus
$\kappa _{{h}} = (k_x^2+k_y^2)^{1/2}$, and the angles
$(\eta _k,\phi _k)$ are such that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn156.png?pub-status=live)
and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn157.png?pub-status=live)
respectively.
The equivalent sources for the particular cases of §§ 4.4 and 5.4 are gathered in table 5. For line sources in two dimensions, obtained as the limit of an ellipse of zero thickness perpendicular to the line, and for plane sources in three dimensions, obtained as the limit of a spheroid of zero thickness perpendicular to the plane, the opposite single layers on either side of the line or plane collapse into a double layer as the thickness goes to zero.
Table 5. Equivalent sources and associated spectra for the oscillating bodies considered in §§ 4.4 and 5.4, with $\varUpsilon$ and
$D(\varUpsilon )$ defined in (5.37)–(5.38).
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_tab5.png?pub-status=live)
For a vibrating circular cylinder, adding viscous attenuation according to (6.23), the same waves are obtained as in Hurley & Keady (Reference Hurley and Keady1997), compared satisfactorily with experiment by Sutherland et al. (Reference Sutherland, Dalziel, Hughes and Linden1999) and Zhang et al. (Reference Zhang, King and Swinney2007). For a vibrating elliptic cylinder and a vertical knife edge, the sources have been shown by Voisin (Reference Voisin2020) to be consistent with the measurements of Sutherland & Linden (Reference Sutherland and Linden2002) and Peacock et al. (Reference Peacock, Echeverri and Balmforth2008), respectively. For a vibrating sphere, the source has been used by Voisin et al. (Reference Voisin, Ermanyuk and Flór2011) to propose an expression of the waves which compares successfully with the experiments of Flynn et al. (Reference Flynn, Onu and Sutherland2003), King et al. (Reference King, Zhang and Swinney2009), Voisin et al. (Reference Voisin, Ermanyuk and Flór2011) and Ghaemsaidi & Peacock (Reference Ghaemsaidi and Peacock2013).
6.3. Relation to the direct approach
Now, the direct approach also leads to an equivalent source (3.13), of spectrum (3.15). Writing, for the elliptic cylinder,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn158.png?pub-status=live)
and using (4.19), we obtain, after some algebra,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn159.png?pub-status=live)
a result significantly more complicated than (6.16). The two spectra coincide at the wave vectors satisfying the dispersion relation $\omega ^2(k_x^2+k_z^2) = N^2k_x^2$, that is
$\sinh \xi _0^2+\sin ^2\eta _k = 0$, namely
$\eta _k = \pm \mathrm {i}\xi _0$, which alone contribute to wave radiation. Accordingly, the effect of the indirect approach is to filter out, in the spectrum of the oscillating body, the wave vectors that do not play a role in its wave radiation.
However, the direct approach is not always possible. This can be seen with the example of an inclined vibrating elliptic cylinder, of major axis making the angle $\phi _0$ to the horizontal. The results of Hurley (Reference Hurley1997), combined with the theory of Voisin (Reference Voisin2020), imply a spectrum of the form given in table 5 for a cylinder of horizontal and vertical axes, but with the angle
$\theta _0$ replaced, in the definition (5.37b) of
$\varUpsilon$, by
$\theta _0+\phi _0$ for the beams of waves propagating upward to the right and downward to the left, and
$\theta _0-\phi _0$ for the beams propagating downward to the right and upward to the left. Such different spectra for the different wave beams cannot be achieved with the indirect approach, and require the switch to the direct approach.
In a related way, Hurley & Keady (Reference Hurley and Keady1997, Appendix B) considered the representation of the cylinder by line distributions of point vortices, and also obtained two different distributions, one per beam, where the vortices are distributed along the critical segment $\zeta _\pm = 0$ joining, for this beam, the critical points on either side of the cylinder; see figures 2 and 3.
The limitation of indirect approaches becomes more explicit for an inclined knife edge, making the angle $\phi _0$ to the horizontal. Introducing inclined coordinates
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn160.png?pub-status=live)
the imposed normal velocity at the knife edge at $z_0 = 0$ induces a pressure discontinuity across it,
$2p_0(x_0) = p(x_0,z_0 = +0)-p(x_0,z_0 = -0)$. Elementary manipulation of the equations of motion for
$0 < \omega < N$ yields a velocity divergence, hence a source of mass,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn161.png?pub-status=live)
When the knife edge is either horizontal ($\phi _0 = 0$) or vertical (
$\phi _0 = {\rm \pi}/2$), a double layer is obtained, consistent with table 5. But when the knife edge is truly inclined, a combination of single and double layers is obtained.
Accordingly, the possibility of using an indirect representation of an oscillating body for generation problems is seen to require some degree of symmetry of the body with respect to the horizontal and the vertical. By contrast, Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) showed that a double layer suffices for scattering problems, irrespective of the shape of the scatterer.
7. Boundary values
The possibility of solving analytically the boundary integral equations for the cylinder and the spheroid has resulted from the availability of expansions (4.18) and (5.15) of the Green's function, allowing the limit to be taken in which $\boldsymbol {x}$ approaches
$S$ asymptotically. When the equations are considered numerically, a procedure is required for which
$\boldsymbol {x}$ is exactly on
$S$. Given the singularity of the Green's function at
$\boldsymbol {x}' = \boldsymbol {x}$, this requires the evaluation of the surface integral as an improper one, in which a vicinity of
$\boldsymbol {x}$ is excluded and the contribution of this vicinity is evaluated separately. Ordinarily, as a result, the Kirchhoff–Helmholtz integral evaluates to one half the value that would be obtained if
$\boldsymbol {x}$ were just outside
$S$; see, for example, Brebbia et al. (Reference Brebbia, Telles and Wrobel1984, chapter 2), Pozrikidis (Reference Pozrikidis1992, chapter 2; Reference Pozrikidis2002, chapters 2 and 4) or Gaul et al. (Reference Gaul, Kögl and Wagner2003, chapter 4).
For monochromatic internal waves, Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) surrounded $\boldsymbol {x}$, as is usually done, by a small hemisphere. They obtained, not the usual factor
$1/2$, but rather a complicated factor depending on the frequency and on the inclination of the normal to
$S$ at
$\boldsymbol {x}$. This is at apparent odds with Gorodtsov & Teodorovich (Reference Gorodtsov and Teodorovich1982), who gave the factor
$1/2$ without justification for the internal waves emitted a body in uniform translation, and Kapitonov (Reference Kapitonov1980), who obtained the same factor
$1/2$ for unsteady inertial waves, based on earlier work by Sobolev (Reference Sobolev1954), replacing the hemisphere by a small vertical cylinder.
7.1. Kirchhoff–Helmholtz integral
We adapt here Sobolev's and Kapitonov's procedure to monochromatic internal waves, on the implicit assumption that $\omega > N$. For this, we revisit first the Kirchhoff–Helmholtz integral (3.4), derived in § 3.1 by distribution theory, applying classical function theory to the volume delimited by the surface
$S$ of the body, the surface
$S_\infty$ at infinity and, if
$\boldsymbol {x}$ is situated outside the body, a circular cylinder
$S_\epsilon$ of centre at
$\boldsymbol {x}$, vertical axis, height
$2h$ and horizontal radius
$\epsilon \to 0$, shown in figure 4(a). For simplicity all lengths are assumed non-dimensional, based on some length scale of the forcing. Sobolev (Reference Sobolev1954) kept
$h$ finite, while Kapitonov (Reference Kapitonov1980) set
$h = \epsilon ^{3/4}$; we keep
$h$ finite for now, and will see later how the necessity arises to set
$h \to 0$ while preserving
$\epsilon /h = o(1)$.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_fig4.png?pub-status=live)
Figure 4. Small cylindrical surface $S_\epsilon$ around a point
$P$ situated either (a) outside or (b) at the surface
$S$, for the calculation of the Kirchhoff–Helmholtz integral at
$P$.
When $\boldsymbol {x}$ is situated inside
$S$, the derivation is exactly the same as before, yielding (3.4b). When
$\boldsymbol {x}$ is situated outside
$S$, we write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn162.png?pub-status=live)
that is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn163.png?pub-status=live)
On $S_\epsilon$ we introduce cylindrical polar coordinates
$(r''_{{h}},\phi '',z'')$ for the position
$\boldsymbol {x}'' = \boldsymbol {x}'-\boldsymbol {x}$ relative to
$\boldsymbol {x}$, as shown in figure 4(a). On the horizontal part
$S_\epsilon ^{({h})}$ of
$S_\epsilon$, we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn164.png?pub-status=live)
yielding
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn165.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn166.png?pub-status=live)
With $\psi$ and
$u_{z''}$ bounded over
$S_\epsilon ^{({h})}$, the contribution of
$S_\epsilon ^{({h})}$ is seen to vanish as
$\epsilon \to 0$. On the vertical part
$S_\epsilon ^{({v})}$ of
$S_\epsilon$, we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn167.png?pub-status=live)
yielding
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn168.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn169.png?pub-status=live)
An additional assumption $h \to 0$ is required, so that not only
$\psi$ and
$u_{r''_{{h}}}$ are bounded on
$S_\epsilon ^{({v})}$, but also
$\psi (\boldsymbol {x}')$ may be made as close to
$\psi (\boldsymbol {x})$ as desired, while ensuring
$\epsilon /h = o(1)$ so that the preceding derivations remain valid. In other words, we require
$\epsilon \ll h \ll 1$, a possible choice being
$h = \epsilon ^{3/4}$, as picked by Kapitonov (Reference Kapitonov1980). The contribution of
$S_\epsilon ^{({v})}$ reduces to
$-\psi (\boldsymbol {x})$, and we finally obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn170.png?pub-status=live)
yielding (3.4a).
When $\boldsymbol {x}$ is situated on
$S$, we write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn171.png?pub-status=live)
where $S_\epsilon$ is the part of the preceding small cylinder situated outside
$S$, up to its intersection with
$S$, delimiting a small domain
$D_\epsilon$ on
$S$, and
$S\setminus D_\epsilon$ is the rest of
$S$ with
$D_\epsilon$ removed, as shown in figure 4(b).
We introduce Cartesian coordinates $(x'',y'',z'')$ on
$S_\epsilon$, with the
$x''$-axis in the vertical plane containing the line of steepest ascent through
$\boldsymbol {x}$ on
$S$, and the
$y''$-axis along the level line. Assuming for definiteness the outward normal
$\boldsymbol {n}$ at
$\boldsymbol {x}$ to be directed upward at the angle
$\theta$ the vertical, the tangent plane has the equation
$z'' = - x''\tan \theta = -r''_{{h}}\tan \theta \cos \phi ''$. The integral over
$S_\epsilon ^{({h})}$ is unchanged, while the integral over
$S_\epsilon ^{({v})}$ comprises two parts: one for
$z''$ varying from
$0$ to
$h$, and the other for
$z''$ varying from
$-r''_{{h}}\tan \theta \cos \phi ''$ to
$0$. This second part cancels out by integration over
$\phi ''$, and we obtain half the integral over the complete small cylinder, namely
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn172.png?pub-status=live)
Alternatively, we might simply point out, as did Kapitonov (Reference Kapitonov1980), that the tangent plane divides the complete cylinder into two parts which are mirror images of each other through the horizontal plane; given the symmetry of the integrand, this implies that the integral over each part is equal to half its value for the complete cylinder.
To evaluate the integral over $S\setminus D_\epsilon$, we adapt the procedure of Martin (Reference Martin2006, § 5.5 and Appendix F). On
$D_\epsilon$, we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn173.png?pub-status=live)
As $\epsilon \to 0$, hence
$r''_{{h}} \to 0$, the domain
$D_\epsilon$ reduces, to leading order, to a small elliptic disc on the tangent plane to
$S$ at
$\boldsymbol {x}$. To the next order, we write
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn174.png?pub-status=live)
so that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn175.png?pub-status=live)
The integral over $D_\epsilon$ is seen to converge at
$r''_{{h}} = 0$ and to vanish as
$\epsilon \to 0$. The integral over
$S\setminus D_\epsilon$ becomes one over the entirety of
$S$ and, though improper, converges at
$\boldsymbol {x}' = \boldsymbol {x}$.
We thus obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn176.png?pub-status=live)
Since the integral is convergent, the factor $1/2$ on the right-hand side is independent of the shape of the small vicinity of
$\boldsymbol {x}$ involved in its derivation, either hemispherical for Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) or cylindrical here. Indeed, the complicated factor derived by Martin & Llewellyn Smith may be shown numerically, if not analytically, to evaluate to
$0.5$ for any values of the frequency ratio
$\omega /N > 1$ and inclination
$\theta$ of the normal
$\boldsymbol {n}$ to the vertical, as seen in Appendix C.
7.2. Layer potentials
This result for the Kirchhoff–Helmholtz integral may be extended to the layer potentials by the method of Courant & Hilbert (Reference Courant and Hilbert1962, § IV.1.4). For this, we consider the double-layer potential (3.10) and continue the dipole density $\mu$ as a twice continuously differentiable function into
$V_+$, subject to the condition that
$(N^2\partial /\partial n_{{h}}-\omega ^2\partial /\partial n)\mu = 0$ on
$S$; we then apply Green's theorem (3.2) inside
$V_+$, choosing
$f = \mu (\boldsymbol {x}')$ and
$g = G(\boldsymbol {x}-\boldsymbol {x}')$. We obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn177.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn178.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn179.png?pub-status=live)
The volume integral on the left-hand side is a continuous function of $\boldsymbol {x}$ across
$S$, implying, as
$\boldsymbol {x}$ approaches
$S$ from either
$V_+$ or
$V_-$, that
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn180.png?pub-status=live)
consistent with (6.6).
We proceed in the same way for the single-layer potential (3.7), continuing the monopole density $\sigma$ as a twice continuously differentiable function into
$V_+$, subject to the condition that
$(N^2\partial /\partial n_{{h}}-\omega ^2\partial /\partial n)\sigma = 0$ on
$S$, then applying (3.2) inside
$V_+$ with
$f = \sigma (\boldsymbol {x}')$ and
$g = G(\boldsymbol {x}-\boldsymbol {x}')$. Taking the derivative
$N^2\partial /\partial n_{{h}}-\omega ^2\partial /\partial n$ of both sides of the result, and letting
$\boldsymbol {x}$ approach
$S$ from either
$V_+$ or
$V_-$, we obtain
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn181.png?pub-status=live)
consistent with (6.3).
8. Far field
At this stage one question remains, the contribution of the surface $S_\infty$ at infinity, assumed to vanish in the derivation of the Kirchhoff–Helmholtz integral (3.4). To answer it, we use the studies of the far field by Voisin (Reference Voisin2003) and Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012).
Consider a sphere $S_R$ of large radius
$R$, with centre within the oscillating body, and evaluate its contribution
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn182.png?pub-status=live)
as $R \to \infty$. The Green's function satisfies
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn183.png?pub-status=live)
In the evanescent frequency range $\omega > N$, the far-field studies give for the layer potentials the decay laws
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn184.png?pub-status=live)
Provided the waves decay at the same rate, and given the surface area element $\mathrm {d}^2S'$ grows as
$R^2$, the contribution of
$S_R$ is seen to vary as
$1/R$ and to vanish as
$R \to \infty$.
In the propagating frequency range $0 < \omega < N$, different decay laws are obtained depending on whether we are inside or outside the wave beams, delimited by the critical rays tangential to
$S$. Outside the beams the decay laws are the same as before, (8.3) for the layer potentials. Provided the waves decay at the same rate, the contribution of that part of
$S_R$ vanishes as
$R \to \infty$. Inside the beams the decay laws become, for the layer potentials,
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn185.png?pub-status=live)
while the associated part of $S_R$ has its area varying as
$R$. Provided the waves decay at the same rate, the contribution of this part of
$S_R$ varies as
$1/R^{1/2}$ and vanishes as
$R \to \infty$.
Accordingly, the decay laws (8.3) and (8.4) must be taken as Sommerfeld-type radiation conditions, to be satisfied by the waves in order to ensure that they all emanate from the oscillating body and do not come in from infinity.
9. Conclusion
Two boundary integral formulations have been presented for the generation of monochromatic internal waves by an oscillating body in an inviscid stratified fluid. One formulation is direct, involves both single and double layers through the Kirchhoff–Helmholtz integral (3.6), and leads to the integral equation (3.5). The other formulation is indirect, involves the single layer (3.7) alone, and leads to the integro-differential equation (6.11). A method for the analytical solution of these equations has been proposed, in two steps: first, in the frequency range of evanescent waves, the horizontal and vertical coordinates are stretched according to (3.17), transforming the wave equation into a Poisson equation to which the usual methods can be applied; and second, the solution is continued analytically to the frequency range of propagating waves, using the causality condition that the waves be analytic, for time dependence as $\exp (-\mathrm {i}\omega t)$ say, in the upper half of the complex
$\omega$-plane, a continuation implemented via Lighthill's radiation condition (3.18).
The method has been applied to the two-dimensional case of an elliptic cylinder in § 4 and the three-dimensional case of a spheroid in § 5. Arbitrary oscillations have been considered, and the results specialized to radial pulsations and rigid vibrations. For these, the results for the waves are gathered in tables 3 and 4, and those for the distributions of singularities equivalent to the body in table 5. The success of the method comes from the existence, in the stretched coordinates, of a separated expansion of the Green's function in eigenfunctions for the given geometry: circular functions for the cylinder, and spherical harmonics for the spheroid. As a result, the boundary integral equation reduces to a diagonal (infinite) linear system, solved analytically. Sturova (Reference Sturova2001), for the cylinder, used a similar expansion but without coordinate stretching; the circular functions were not eigenfunctions, yielding a non-diagonal system solved numerically by truncation.
The waves for arbitrary oscillations of the cylinder are obtained as an expansion of the general form anticipated by Barcilon & Bleistein (Reference Barcilon and Bleistein1969) and Hurley (Reference Hurley1972). They involve complex stretched coordinates, which Lighthill's radiation condition allows to be written as simple combinations of the characteristic coordinates $x_+$ and
$x_-$, with the determination of these multivalued combinations set unequivocally.
The greatest interest of the boundary integral method lies, however, elsewhere; namely, in the determination of a distribution of singularities equivalent to the oscillating body. First, this distribution allows the effect of viscosity to be added into the theory. Inviscid waves are but a web of singular critical rays, at which the amplitude diverges and the phase jumps; it is viscosity, acting locally to smooth out the singularities, which determines the distribution of the waves into space and sets ultimately their amplitude and their phase. The mathematical description of this process takes place in Fourier space, as shown by Lighthill (Reference Lighthill1978, § 4.10) at large distances from the forcing and Voisin (Reference Voisin2020) at arbitrary distance from it. The knowledge of a representation of the body as a distribution of singularities, and thence of its spectrum, provides the basis for this description; it avoids the need, as did Hurley (Reference Hurley1997) and Hurley & Keady (Reference Hurley and Keady1997), to calculate the inviscid waves first, then derive their Fourier representation and finally add viscous attenuation into it.
The representations of the vibrating elliptic cylinder and sphere in table 5 have been applied by Voisin (Reference Voisin2020) and Voisin et al. (Reference Voisin, Ermanyuk and Flór2011), respectively, to the calculation of their waves and compared successfully with experiment. The effect of viscosity is illustrated in figure 5 for the vertical vibrations at the frequency $\omega$ and velocity
$W = -\omega A$, with
$\omega _0/N = 0.44$ and
$A = 0.32\ \mathrm {cm}$, of a cylinder of horizontal semi-axis
$a = 2.52\ \mathrm {cm}$, vertical semi-axis
$b = 0.86\ \mathrm {cm}$ and aspect ratio
$a/b \approx 3$ in a fluid of buoyancy frequency
$N = 0.97\ \mathrm {s}^{-1}$ and kinematic viscosity
$\nu = 1\ \mathrm {mm}^2\ \mathrm {s}^{-1}$, corresponding to the measurements in figure 8(d) of Sutherland & Linden (Reference Sutherland and Linden2002). The plotted quantity is the amplitude of the time derivative of the buoyancy frequency disturbance,
$N_t^2 = -N^2\partial u_z/\partial z$, normalized by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn186.png?pub-status=live)
The inviscid expression for this quantity, derived from table 3, is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn187.png?pub-status=live)
and its viscous expression, derived from (6.23), is
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn188.png?pub-status=live)
Both expressions are plotted in figure 5. To better illustrate the structure of the waves, in particular the occurrence of singularities at critical points on the cylinder, they are plotted not only in the fluid but also inside the body. Inviscid dynamics is seen to provide a wave skeleton, which is then fleshed out by viscous diffusion.
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_fig5.png?pub-status=live)
Figure 5. Waves for the elliptic cylinder in figure 8(d) of Sutherland & Linden (Reference Sutherland and Linden2002), as predicted (a) in an inviscid fluid by (9.2) and (b) in a viscous fluid by (9.3). The predictions are plotted both outside and inside the cylinder, whose outline is shown dashed.
Another application of the distribution of singularities is, for a rigid body, the determination of the added mass of the body. The added mass follows from the first moment of the distribution, and provides direct access to the energy radiated by the body and to the force exerted on it, without requiring the actual calculation of the waves. The energy is important for oceanic applications, as it characterizes the conversion rate from the barotropic tide oscillating over bottom topography to the internal tide that this oscillation generates; and the force is important for the stability of midwater floats, as it governs their oscillation back to their neutral buoyancy level once displaced away from it then released. These aspects will be considered separately; a summary presentation may be found in Voisin (Reference Voisin2009).
Coming back to the foundations of the boundary integral method, it must be noted that the derivations in §§ 3 and 7 have all been made for evanescent waves, of frequency $\omega > N$. For these, the only singularities in the boundary integrals arise when the observation point
$P$ is on the surface
$S$. For propagating waves, of frequency
$0 < \omega < N$, the situation is different: singularities are present whenever
$P$ is inside a wave beam, delimited by the critical rays tangential to
$S$, and they take place at one or several closed curves where the double cone with apex at
$P$, vertical axis and semi-aperture
$\theta _0$ intersects
$S$. Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) considered this frequency range in detail, paying particular attention to the far field, to which the present § 8 is devoted. It is beyond our scope to investigate this frequency range further. As a rule, we expect those of our conclusions that do not involve the frequency explicitly to hold whatever the frequency.
One such conclusion concerns the possibility of replacing the direct approach of § 3, representing the oscillating body as a combination of single and double layers, by the indirect approach of § 6, involving a single layer alone. For generation problems, this possibility requires the existence of a fictitious flow inside the body, imagined to be filled with fluid, simultaneously satisfying the wave equation throughout the body and the continuity of the pressure at its boundary. The elliptic cylinder and the spheroid in §§ 4 and 5, symmetric with respect to the horizontal and vertical directions, may be treated by both approaches; an elliptic cylinder with inclined axes requires, however, the direct approach. For scattering problems, Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) showed that a double layer always suffices.
Another conclusion concerns the values of the boundary integrals at the boundary. For the Kirchhoff–Helmholtz integral (7.13), the same arithmetic mean $\psi (\boldsymbol {x})/2$ between the values
$\psi (\boldsymbol {x})$ and
$0$ of the integral (3.4) on either side of the boundary is obtained as for the Laplace and Helmholtz equations. Similarly, the same additional term
$\pm \mu (\boldsymbol {x})/2$, with
$\mu (\boldsymbol {x})$ the dipole density, is obtained for the double-layer potential (7.15) when the boundary is approached from either side, and the same additional term
$\pm \sigma (\boldsymbol {x})/2$, with
$\sigma (\boldsymbol {x})$ the monopole density, for the modified normal derivative (7.16) of the single-layer potential.
It is unclear at this stage how boundary elements would be implemented numerically for internal waves in a three-dimensional setting. All existing numerical implementations are two-dimensional and deal with vertically trapped waves, decomposed into a series of normal modes. Sturova (Reference Sturova2006, Reference Sturova2011) considers a cylinder, either circular or elliptic, and uses an indirect single-layer representation for $0 < \omega < N$ and a direct Kirchhoff–Helmholtz representation for
$\omega > N$; she does not discuss, however, how the contour of the cylinder is decomposed into elements, nor the way the weak singularities of the kernels are processed. Pétrélis et al. (Reference Pétrélis, Llewellyn Smith and Young2006) consider various topographies and use straight elements, discretized with respect to the vertical coordinate; Balmforth & Peacock (Reference Balmforth and Peacock2009), Echeverri & Peacock (Reference Echeverri and Peacock2010), Echeverri et al. (Reference Echeverri, Yokossi, Balmforth and Peacock2011) and Mathur et al. (Reference Mathur, Carter and Peacock2014) do the same but discretize the elements with respect to the horizontal coordinate. The topographies are represented as single layers and the weak singularities of the kernels are eliminated by integrating the integral equations along each element. Clearly, the boundary element method for three-dimensional internal waves remains a topic for further investigation.
Acknowledgements
This research was initiated at DAMTP under the supervision of D. Crighton in 1998, and continued at LEGI where it benefitted from the input of A. Chikhi during a master internship in 2005 and from conversations with S. Llewellyn Smith and A. Davis during a stay at UCSD in 2008. N. Shmakova is thanked for kindly translating the obituary by Sveshnikov et al. (Reference Sveshnikov, Shishmarev and Pletner1989), and E. Ermanyuk for clarifying the translation of Sobolev (Reference Sobolev1954). Application of the boundary integral method to monochromatic internal waves was suggested by V. Gorodtsov at the Institute for Problems in Mechanics in Moscow on 26 May 1998, during a marathon conversation organized by Y. Chashechkin with S. Smirnov acting as a translator. One anonymous reviewer is thanked for pointing out several misconceptions about boundary integrals in the original manuscript and, together with the editor, for thoughtful suggestions which led to redesign and improvement of the paper.
Funding
This work was supported by Marie Curie grant ERBFMBI CT97 2653 (1998) from the European Commission, grants TOPOGI-3D (2005–2008) and PIWO (2008–2011) from the French National Research Agency (ANR), and a grant from Fondation Simone et Cino Del Duca–Institut de France (2015–2018).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Legendre functions
The ordinary Legendre functions $P_l(z)$ and
$Q_l(z)$ of non-negative integer degree
$l$ are defined in the plane of the complex variable
$z$, cut along the real interval
$[-1,1]$, as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn189.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn190.png?pub-status=live)
The associated Legendre functions $P_l^m(z)$ and
$Q_l^m(z)$ are defined for non-negative integer order
$m$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn191.png?pub-status=live)
and for negative order as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn192.png?pub-status=live)
On the cut $[-1,1]$, the Legendre functions are also called Ferrers functions. The ordinary functions
$\operatorname {P}_l(x)$ and
$\operatorname {Q}_l(x)$ of the real variable
$x \in [-1,1]$ are defined for degree
$l$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn193.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn194.png?pub-status=live)
and the associated functions $\operatorname {P}_l^m(x)$ and
$\operatorname {Q}_l^m(x)$ are defined for positive order
$m$ as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn195.png?pub-status=live)
and for negative order as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn196.png?pub-status=live)
The two types of functions are related to each other at the cut by
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn197.png?pub-status=live)
Several definitions and notations of the Legendre functions exist in the literature. The above definitions and notations are based on Hobson (Reference Hobson1931, chapters 2, 3 and 5), Abramowitz & Stegun (Reference Abramowitz and Stegun1972, chapter 8) and Olver et al. (Reference Olver, Lozier, Boisvert and Clark2010, chapter 14).
The derivations in § 5 use the symmetry relations
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn198.png?pub-status=live)
and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn199.png?pub-status=live)
the recurrence relations
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn200.png?pub-status=live)
and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn201.png?pub-status=live)
with $Q_l^m(z)$ satisfying the same relations as
$P_l^m(z)$, and
$\operatorname {Q}_l^m(x)$ as
$\operatorname {P}_l^m(x)$, and the identities
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn202.png?pub-status=live)
and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn203.png?pub-status=live)
Gegenbauer's finite integral in § 12.14 of Watson (Reference Watson1944) is also used, in the form
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn204.png?pub-status=live)
with $\operatorname {J}_m(z)$ a cylindrical Bessel function and
$\operatorname {j}_l(z)$ a spherical Bessel function.
Appendix B. Spherical harmonics
Spherical harmonics are defined in terms of associated Legendre functions as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn205.png?pub-status=live)
with $l$ a non-negative integer and
$m$ an integer such that
$|m| \leqslant l$; see for example Jackson (Reference Jackson1999, § 3.5) and Olver et al. (Reference Olver, Lozier, Boisvert and Clark2010, § 14.30). They form a complete orthonormal set over the surface of the unit sphere, in terms of which an arbitrary function
$g(\theta ,\phi )$ may be expanded as
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn206.png?pub-status=live)
with coefficients
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn207.png?pub-status=live)
where an overbar denotes a complex conjugate and $\mathrm {d}\varOmega = \sin \theta \,\mathrm {d}\theta \,\mathrm {d}\phi$ the solid angle element.
Appendix C. Another look at boundary values
When $\boldsymbol {x}$ is on
$S$ and the Kirchhoff–Helmholtz integral is evaluated on a small surface surrounding it, Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) used a hemisphere
$H_\epsilon$ of radius
$\epsilon$ with its polar axis along the normal
$\boldsymbol {n}$ to
$S$ at
$\boldsymbol {x}$. Introducing spherical polar coordinates
$(r'',\theta '',\phi '')$ with the same polar axis (so that
$\phi ''$ differs from that in § 7.1), we have
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn208.png?pub-status=live)
yielding
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn209.png?pub-status=live)
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn210.png?pub-status=live)
where $\varGamma = N/\omega$, with
$0 < \varGamma < 1$, and
![](https://static.cambridge.org/binary/version/id/urn:cambridge.org:id:binary:20210920093657605-0327:S0022112021007291:S0022112021007291_eqn211.png?pub-status=live)
As $\epsilon \to 0$ the same result (7.9) is obtained as for the cutoff vertical cylinder
$S_\epsilon$ in figure 4(b), but with the factor
$1/2$ on the right-hand side replaced by the right-hand side of (C2b).
Now, as seen in § 7.1, the two integrals over $S_\epsilon$ and
$H_\epsilon$ must have the same limit as
$\epsilon \to 0$. Numerically, the right-hand side of (C2b) evaluates to
$0.5$ whatever the values of
$\varGamma$ and
$\theta$. Analytically, it evaluates to
$1/2$ for
$\varGamma = 0$, independent of
$\theta$, but the extension of this result to arbitrary
$\varGamma$ and
$\theta$ has not been possible. We may, however, remark, as for
$S_\epsilon$ following (7.9), that the tangent plane to
$S$ at
$\boldsymbol {x}$ divides a sphere of radius
$\epsilon$ around
$\boldsymbol {x}$ into two hemispheres which are mirror images of each other through the horizontal plane, and that the symmetry of the integrand is such that each hemisphere gives half the integral over the complete sphere, shown by Martin & Llewellyn Smith (Reference Martin and Llewellyn Smith2012) to evaluate to
$1$.