Spheroid and Disc
The oblate and the prolate spheroid, as well as the disc, a degenerate oblate spheroid, are solved by a series in the spheroidal wave functions of their shape, which is shared by every acoustic excitation. Both the plane wave and the monopole are supported, at an arbitrary direction of incidence or position, and the scatterer may have an arbitrary orientation. How the solution is obtained, and how accurate it is, is described under Spheroidal Solution.
Geometry
The scatterer is a surface $\xi = \xi_0$ of spheroidal coordinates $(\xi, \eta, \varphi)$. It is specified by its equatorial radius $a$ and its polar radius $b$, the latter measured along the axis of revolution, and the shape follows from them:
| radii | semifocal distance $f$ | surface | coordinates | |
|---|---|---|---|---|
OblateSpheroid | $a > b$ | $\sqrt{a^2 - b^2}$ | $\xi_0 = b / f \ge 0$ | oblate |
ProlateSpheroid | $b > a$ | $\sqrt{b^2 - a^2}$ | $\xi_0 = b / f > 1$ | prolate |
The special case $b = 0$ of the oblate spheroid, that is, $\xi_0 = 0$, is the disc of radius $a = f$, whose two faces are distinguished by the sign of $\eta$. The corresponding degeneration of the prolate spheroid, $a = 0$ or $\xi_0 = 1$, is the segment of the axis between the foci, which is not supported.
$a \neq b$ is required: the spheroidal coordinates degenerate for a sphere, where $f \rightarrow 0$ and $\xi_0 \rightarrow \infty$. A sphere is a HardSphere or a SoftSphere, for which the considerably cheaper spherical series is evaluated.
The boundary condition is carried as a type parameter, as it is orthogonal to the geometry:
| type | condition on the surface |
|---|---|
SoundHard | $\gamma_1 p = \hat{\bm n} \cdot \nabla p = 0$ (the normal velocity vanishes) |
SoundSoft | $\gamma_0 p = 0$ (pressure release) |
API
SphericalScattering.Spheroid — Type
Spheroid{BC}(
equatorialRadius = error("missing argument `equatorialRadius`"),
polarRadius = error("missing argument `polarRadius`"),
axis = SVector(0.0, 0.0, 1.0)
)Constructor for a spheroid centered in the origin, where BC is SoundHard or SoundSoft.
The spheroid is specified by its equatorialRadius $a$ and its polarRadius $b$, the latter being measured along the axis of revolution, which is normalized. The shape follows from the radii: for $a > b$ an OblateSpheroid is returned, with the semifocal distance $f = \sqrt{a^2 - b^2}$ and the radial coordinate $ξ_0 = b / f$ of the surface; for $b > a$ a ProlateSpheroid, with $f = \sqrt{b^2 - a^2}$ and again $ξ_0 = b / f$.
For polarRadius = 0 the oblate spheroid degenerates into a disc of radius $a$, see Disc. The prolate spheroid degenerating into a line segment, equatorialRadius = 0, is not supported.
$a ≠ b$ is required: the spheroidal coordinates degenerate for a sphere, as $f → 0$ and $ξ_0 → ∞$ in that limit. Use HardSphere or SoftSphere for a sphere.
SphericalScattering.Disc — Function
Disc(BC; radius = error("missing argument `radius`"))Constructor for a disc of radius $a$, that is, for the degenerate oblate spheroid with vanishing polar radius, where BC is SoundHard or SoundSoft.
Its two faces are distinguished by the sign of $η$; the surface is $ξ_0 = 0$.
The axis of revolution $\hat{\bm a}$ is normalized during the initialization and defaults to $\hat{\bm e}_z$. Since the pressure is a scalar, an arbitrarily oriented spheroid costs no more than transforming the evaluation points into the frame of the scatterer; nothing has to be rotated back.
The dimensions can be queried, which is of interest mainly because a Spheroid stores $f$ and $\xi_0$ rather than the radii it was constructed from:
SphericalScattering.equatorialRadius — Function
equatorialRadius(sp::Spheroid)Returns the equatorial radius of the spheroid, $f \sqrt{1 + ξ_0^2}$ for an oblate and $f \sqrt{ξ_0^2 - 1}$ for a prolate one.
SphericalScattering.polarRadius — Function
polarRadius(sp::Spheroid)Returns the polar radius $f ξ_0$ of the spheroid, which vanishes for a disc.
SphericalScattering.isdisc — Function
isdisc(sp::Spheroid)Returns whether the spheroid is degenerate, that is, a disc.
Quantities
The general API is employed, with or without the modal coefficients:
p = scatteredfield(sp, ex, Pressure(point_cart))
FF = scatteredfield(sp, ex, FarField(point_cart))
γ₀ = scatteredfield(sp, ex, PressureTrace(point_cart))
γ₁ = scatteredfield(sp, ex, PressureNormalGradient(sp, point_cart))
Δp = scatteredfield(sp, ex, PressureJump(point_cart)) # a disc onlyand likewise field for the total quantities, where sp is a Spheroid — possibly the degenerate one returned by Disc. As for a sphere, the total far field exists for a monopole but not for a plane wave, and the pressure is set to zero inside the scatterer, where the oblate spheroidal coordinates answer what the radius of a sphere cannot.
The far field omits the factor $\mathrm{e}^{-\mathrm{j} k r} / r$, as for the other scatterers, see Quantities; how the quantities follow from the modes is described under The Quantities in Terms of the Modes.
Surface Traces
The traces are evaluated on the surface: only the direction of each location is taken into account, its radial coordinate being replaced by $\xi_0$. Hence the locations may also be given by the points of a faceted surface mesh, which do not lie exactly on the spheroid.
The outward normal of a spheroid is not parallel to the position vector — for a disc it is $\pm \hat{\bm e}_z$ everywhere. The default of PressureNormalGradient, the radial direction, is therefore wrong for a spheroid. Pass the scatterer, which fills in the correct normals:
γ₁ = scatteredfield(sp, ex, PressureNormalGradient(sp, point_cart))SphericalScattering.outwardNormal — Function
outwardNormal(sphere::Spheroid, point)Compute the outward unit normal of the spheroid at the surface point having the direction of point.
The normal is $\hat{e}_ξ$ of the spheroidal coordinates, which is not parallel to the position vector unless the scatterer is a sphere: for a disc it is $\pm \hat{e}_z$, the sign following the face. The point and the returned normal are in Cartesian coordinates of the global frame.
SphericalScattering.outwardNormals — Function
outwardNormals(sphere::Spheroid, locations)Compute the outward unit normals of the spheroid at all locations, see outwardNormal.
The Jump Across a Disc
On an open surface a one-sided trace is not determined by the geometry alone, which is why the natural unknown of a boundary element formulation is the jump $[p] = p|_+ - p|_-$, see PressureJump. The two faces of the degenerate surface $\xi = 0$ are $\eta > 0$ and $\eta < 0$. The jump of the total field equals that of the scattered field, the incident field being continuous across the disc.
At the degenerate surface the radial functions split by parity, and the two boundary conditions select opposite halves of the spectrum: the Neumann problem is carried by the modes with odd $n - m$, the Dirichlet problem by those with even $n - m$. Consequently $[p]$ is carried entirely by the sound-hard disc and vanishes identically for the sound-soft one, whose unknown is the jump of the normal derivative instead — that quantity is not implemented yet.
The same parity governs the one-sided traces, for which the face with $\eta > 0$ is returned:
| disc | $\gamma_0$ | $\gamma_1$ |
|---|---|---|
| sound-soft | same on both faces | opposite |
| sound-hard | opposite | vanishes |
Reusing the Modal Coefficients
The spheroidal wave functions are expensive compared to the spherical ones. The coefficients are therefore computed once for a given scatterer and excitation and are reused for every evaluation point. They are also accessible directly, in order to be reused across several calls:
md = SphericalScattering.modes(sp, ex) # the expensive part
p = scatteredfield(sp, ex, md, Pressure(points)) # cheap from here on
FF = scatteredfield(sp, ex, md, FarField(directions))SphericalScattering.modes — Function
modes(sphere::Spheroid, excitation::AcousticExcitation; M=nothing, N=nothing, parameter::Parameter=Parameter())Compute the modal coefficients of the scattering problem, see SpheroidalModes.
The incident field is expanded as $p_\mathrm{i} = \sum_{mn} A_{mn} R^{(1)}_{mn}(c, ξ) S_{mn}(c, η) \mathrm{e}^{\mathrm{j}mφ}$, the coefficients $A_{mn}$ being known analytically for the plane wave and the monopole, see incidentCoefficients. Since the surface $ξ = ξ_0$ is a coordinate surface, the boundary condition decouples the modes, each being scattered with the coefficient $b_{mn}$ of scatterCoeff.
The truncations are determined automatically unless they are given. Since the scattering coefficients decay once a mode is cut off at the surface, the degree is bounded by the size of the scatterer, so that $N = \lceil c \, ρ \rceil + 15$ is taken initially, where $ρ$ is the radius of the circumscribing sphere in units of the semifocal distance, see normalizedCircumradius; $c ρ$ is the counterpart of $ka$. The order M is not estimated but measured, see modesForDegree: a nearly axial excitation needs far fewer orders than a grazing one, which matters, as every order costs calls to the spheroidal wave functions.
The truncation is verified rather than trusted, on the surface of the scatterer: the size of a mode of the scattered pressure there bounds its contribution everywhere, as the outgoing radial functions decrease outward. The modes of the last degrees have to be negligible. A source close to the scatterer requires more degrees than its size suggests, as the expansion of its field converges the more slowly the closer the source is. If the degree has been determined automatically, it is therefore raised by half until the omitted modes are negligible, up to four times its initial value and as long as the result keeps improving: for a source very close to the scatterer the radial functions of high degree eventually exceed the range of the floating-point numbers, which limits the attainable accuracy. A message is printed if the relative accuracy is not attained. A given N, or a nmax of Parameter, is used as it is.
What is stored is the spheroidal parameter $c$, the two truncations, and the two coefficient tables $A_{mn}$ and $b_{mn}$, indexed as [m + M + 1, n + 1], see Solution Approach for their meaning, and see SpheroidalModes.
The convenience signatures without md determine the coefficients internally and discard them afterwards. Whenever more than one quantity is of interest for the same scatterer and excitation, passing md explicitly avoids recomputing them.
Limitations
- A sphere is not a limiting case that can be evaluated: the coordinates degenerate, hence $a \neq b$ is enforced. Use
HardSphereorSoftSphere. - The prolate spheroid degenerating into a line segment, $a = 0$, is not supported; a thin one close to it is, and is as accurate as any other.
- The jump of the normal derivative across a disc, the natural unknown of the sound-soft case, is not implemented.
- A monopole has to lie outside the scatterer, as the expansion of its field assumes. This is checked, an error being thrown otherwise.