1 Introduction
This article studies planar self-similar incompressible vortex sheets meeting each other or a solid wall in cusps (figure 3 b, figure 1 b, figure 8 a,b). Here self-similar flows, also called ‘pseudo-steady’ or ‘quasi-steady’, have vorticity
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D714}(t,\boldsymbol{x})=t^{-1}\unicode[STIX]{x1D714}(1,t^{-\unicode[STIX]{x1D707}}\boldsymbol{x}) & & \displaystyle\end{eqnarray}$$
for some similarity exponent
$\unicode[STIX]{x1D707}$
. The flow is ‘similar’ at all times
$t>0$
, with any spatial distance
$L$
dilated by a factor
$t^{\unicode[STIX]{x1D707}}$
and any velocity
$\boldsymbol{v}$
scaled accordingly. Such flows arise from initial data of type
$$\begin{eqnarray}\unicode[STIX]{x1D714}(0,\boldsymbol{x})=r^{-\unicode[STIX]{x1D705}}{\unicode[STIX]{x1D714}\unicode[STIX]{x0030A}}(\unicode[STIX]{x1D703}),\end{eqnarray}$$
where
$r,\unicode[STIX]{x1D703}$
are polar coordinates. Since vorticity has dimensions of inverse time, the exponent
$-\unicode[STIX]{x1D705}$
of the initial singularity fixes
$t^{-1}\sim L^{-\unicode[STIX]{x1D705}}$
and thus determines the similarity exponent
$\unicode[STIX]{x1D707}=1/\unicode[STIX]{x1D705}$
. To ensure velocity is locally integrable at the origin we only consider
$\unicode[STIX]{x1D707}>{\textstyle \frac{1}{2}}$
. The case
$\unicode[STIX]{x1D707}=1$
is commonly studied for compressible flow as well since it corresponds to
$$\begin{eqnarray}\displaystyle \boldsymbol{v}(t,\boldsymbol{x})=\boldsymbol{v}(1,t^{-1}\boldsymbol{x}). & & \displaystyle\end{eqnarray}$$
If
$\unicode[STIX]{x1D714}$
is single signed, or if the other sign is present but small enough to be ‘dominated’, then strong rotation at the origin winds
$\unicode[STIX]{x1D714}$
(in particular its
$\unicode[STIX]{x1D714}=0$
level sets) into algebraic spirals of type
$r\sim \unicode[STIX]{x1D703}^{-\unicode[STIX]{x1D707}}$
(figure 2
b). For vortex sheets algebraic spiral flows have been modelled in Kaden (Reference Kaden1931), Rott (Reference Rott1956), Stern (Reference Stern1956), Birkhoff (Reference Birkhoff1962), Mangler & Weber (Reference Mangler and Weber1967) and Moore (Reference Moore1975); recently, Elling (Reference Elling2016) gave a mathematically rigorous existence proof for a class of mixed-sign smooth
${\unicode[STIX]{x1D714}\unicode[STIX]{x0030A}}$
. (For logarithmic spirals with entirely different scaling see e.g. Prandtl (Reference Prandtl1924), Alexander (Reference Alexander1971), Kaneda (Reference Kaneda1989) and Elling & Gnann (Reference Elling and Gnann2019).)
It is natural to ask what happens in the borderline case where neither sign of vorticity dominates. In this article a pair of vortex sheets of equal strength but opposite sign is considered. It is shown that in some cases they may form a spiral-free cusp of type
$y\sim \unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FC}}$
(with
$y$
the distance to the symmetry axis,
$\unicode[STIX]{x1D709}$
the distance to the cusp along the axis); otherwise a spiral-ended jet appears.
Figure 1. (a) Initial data with
$N=4$
symmetric pairs of vortex sheets (blue/solid
$\unicode[STIX]{x1D714}>0$
, red/dashed
$\unicode[STIX]{x1D714}<0$
); (b) for
$t>0$
each pair forms a cusp with cusp point at the origin.
Figure 2. (a) Initial data with mixed-sign
$\unicode[STIX]{x1D714}$
, contours are
$\unicode[STIX]{x1D714}=0$
level sets separating
$\unicode[STIX]{x1D714}\gtrless 0$
sectors. (b) Flow for
$t>0$
; (c)
$t\gg 0$
.
Figure 3. Self-similar simple Mach reflection of a shock at a solid corner producing a vortex sheet. Jetting (c) is favoured for shallow ramps and large incident shock Mach numbers, otherwise the sheet appears to blend into the wall in a cusp (b); see Van Dyke (Reference Van Dyke1982, figure 325 etc.) for experiments.
While ‘equal strength’ may at first glance appear to be a less important borderline case, it is precisely the case of a self-similar vortex sheet merging into a slip condition wall (figure 3 b), since eliminating the wall by a reflection yields a mirror sheet of opposite vorticity. Vortex–wall merging is commonly observed, for example in Mach reflection (figure 3; see Mach & Wosyka (Reference Mach and Wosyka1875), von Neumann (Reference von Neumann1943), Hornung (Reference Hornung1986) and Ben-Dor (Reference Ben-Dor1992); experimental observations e.g. Van Dyke (Reference Van Dyke1982) figures 236, 237). Henderson et al. (Reference Henderson, Vasilev, Ben-Dor and Elperin2003) and Vasilev et al. (Reference Vasilev, Ben-Dor, Elperin and Henderson2004) have begun a theoretical and numerical study on when the sheet merges smoothly rather than flipping backwards and forming a spiral-ended jet near the wall.
In § 2 we discuss several heuristic ordinary differential equation (ODE) models derived from the Birkhoff–Rott equation. Section 2.2 approximates the axis-parallel on-sheet velocity
$v^{x}$
as a function of
$\unicode[STIX]{x1D709}$
alone, deriving a formula
$v^{x}\sim \unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FD}-1}$
(see (2.38)) with exponent
$\unicode[STIX]{x1D6FD}$
depending only on
$\unicode[STIX]{x1D707}$
and strain parameter
$e_{1}$
, a single real number representing the limit of
$\unicode[STIX]{x2202}v^{x}/\unicode[STIX]{x2202}\unicode[STIX]{x1D709}$
outside the cusp. Section 2.3 derives a model for the cusp shape as
$y\sim \unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FC}}$
, predicting cusp exponents
$\unicode[STIX]{x1D6FC}$
which are also functions of
$\unicode[STIX]{x1D707},e_{1}$
alone (see (2.45)).
In § 3 we test our model against numerical data. Section 3.1 describes the numerics used for calculating cusps and the new obstacles compared to previous sheet calculations. Section 3.2 compares piecewise linear reconstruction to quadratic and higher degree; while the latter produce matching results that agree with our preferred model, we find that linear reconstruction produces a cusp exponent matching oversimplified conservation-violating models. (This is in contrast to the case of vortex spirals for which even point vortex methods were successful (Pullin Reference Pullin1978, Reference Pullin and Caflisch1989).) Sections 3.2.3 and 3.2.4 confirm that our
$\unicode[STIX]{x1D6FD},\unicode[STIX]{x1D6FC}$
formulas agree with numerical observations, with correct dependence on
$\unicode[STIX]{x1D707}$
and
$e_{1}$
. Section 3.2.5 shows that clockwise circulation (on the upper sheet for a cusp opening to the right) is sometimes possible if
$e_{1}$
is sufficiently large, although generally counterclockwise circulation is observed, which matters because that is the circulation produced by Mach reflections (see von Neumann Reference von Neumann1943; Henderson & Menikoff Reference Henderson and Menikoff1998; Serre Reference Serre2007; Elling Reference Elling2019a
). In cases where cusps do not appear we generally observe spiral-ended jets instead, as in § 3.2.6 where we find that, even for parameters far from the jet–cusp transition, the jets can be rather small, mimicking cusps, suggesting that jetting is frequently obscured by viscous or kinetic effects in numerics or experiments, or simply small enough to be overlooked.
We conclude that our model is correct, that vortex cusps exist and can be computed accurately using vortex segments of at least quadratic degree.
2 Modelling
2.1 Birkhoff–Rott equation
Figure 4. (a) A point
$z^{\prime }$
on the sheet and its mirror image
${z^{\prime }}^{\ast }$
induce almost opposite velocities in a distant point
$z_{0}$
. (b) Two points near
$z_{0}$
and at equal distances induce almost opposite velocities.
We are interested in a vortex sheet in the upper half-plane meeting the wall in a ‘cusp’. We may assume the sheet approaches the cusp point from the right since approach from the left reduces to this case by mirror reflection. Equivalently, we may think of a pair of vortex sheets in the whole plane (figure 4), mirror symmetric across the horizontal axis, one above one below, with circulation of opposite sign. Here,
$\tilde{\unicode[STIX]{x1D6E4}}$
is the circulation on the upper sheet,
$\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}})$
with
$\tilde{z}=\tilde{x}+\text{i}{\tilde{y}}$
the location of the corresponding sheet point at time
$t$
. The infinitesimal segment
$\tilde{\unicode[STIX]{x1D6E4}}^{\prime }$
to
$\tilde{\unicode[STIX]{x1D6E4}}^{\prime }+\text{d}\tilde{\unicode[STIX]{x1D6E4}}^{\prime }$
of the upper sheet, regarded as an infinitesimal point vortex, induces at another point
$\tilde{z}_{0}$
the complex velocity
$$\begin{eqnarray}\displaystyle \tilde{v}^{x}-\text{i}\tilde{v}^{y}=\frac{1}{2\unicode[STIX]{x03C0}\text{i}}\frac{\text{d}\tilde{\unicode[STIX]{x1D6E4}}^{\prime }}{\tilde{z}_{0}-\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}}^{\prime })}. & & \displaystyle\end{eqnarray}$$
The entire upper sheet induces a complex velocity
$$\begin{eqnarray}\displaystyle w(\tilde{z}_{0})=\text{p.v.}\int \frac{\text{d}\tilde{\unicode[STIX]{x1D6E4}}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}(\tilde{z}_{0}-\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}}^{\prime }))}, & & \displaystyle\end{eqnarray}$$
where ‘p.v.’ indicates principal value, the arithmetic average of the velocity limits on each side, i.e. of
$w(\tilde{z}_{\pm })$
as
$\tilde{z}_{\pm }$
approach
$\tilde{z}_{0}$
from each side of the sheet. For the lower sheet, contribution
$\text{d}\tilde{\unicode[STIX]{x1D6E4}}^{\prime }$
is replaced with its negative (opposite circulation),
$\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}}^{\prime })$
with its conjugate (mirror image across the horizontal axis). The upper sheet evolves according to the Birkhoff–Rott (BR) equation
$$\begin{eqnarray}\displaystyle \tilde{z}_{t}(t,\tilde{\unicode[STIX]{x1D6E4}})=\tilde{w}^{\ast },\quad \tilde{w}=\text{p.v.}\int \frac{1}{\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}})-\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}}^{\prime })}-\frac{1}{\tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}})-\tilde{z}^{\ast }(t,\tilde{\unicode[STIX]{x1D6E4}}^{\prime })}\frac{\text{d}\tilde{\unicode[STIX]{x1D6E4}}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}},\quad & & \displaystyle\end{eqnarray}$$
where subscripts indicate partial derivatives. The integration domain is chosen based on the case at hand. The BR equation is equivalent to the incompressible Euler equations under mild assumptions (see e.g. Lopes-Filho, Lopes & Schochet (Reference Lopes-Filho, Lopes and Schochet2007)).
We are interested in self-similar vortex sheets: since circulation
$\tilde{\unicode[STIX]{x1D6E4}}$
has dimensions of length squared over time, it is scaled by a factor
$t^{2\unicode[STIX]{x1D707}-1}$
. We pass to dimensionless
$z,\unicode[STIX]{x1D6E4}$
with the ansatz
$$\begin{eqnarray}\displaystyle \tilde{z}(t,\tilde{\unicode[STIX]{x1D6E4}})=t^{\unicode[STIX]{x1D707}}z(1,t^{1-2\unicode[STIX]{x1D707}}\unicode[STIX]{x1D6E4}), & & \displaystyle\end{eqnarray}$$
resulting in the self-similar Birkhoff–Rott equation
$$\begin{eqnarray}0=w^{\ast }-\unicode[STIX]{x1D707}z+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}z_{\unicode[STIX]{x1D6E4}},\quad w=\text{p.v.}\int \frac{1}{z(\unicode[STIX]{x1D6E4})-z(\unicode[STIX]{x1D6E4}^{\prime })}-\frac{1}{z(\unicode[STIX]{x1D6E4})-z^{\ast }(\unicode[STIX]{x1D6E4}^{\prime })}\frac{\text{d}\unicode[STIX]{x1D6E4}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}}.\end{eqnarray}$$
Here,
$w^{\ast }-\unicode[STIX]{x1D707}z$
is the complex form of the pseudo-velocity
$\boldsymbol{q}=\boldsymbol{v}-\unicode[STIX]{x1D707}\boldsymbol{x}$
;
$\unicode[STIX]{x1D6E4}\mapsto z(\unicode[STIX]{x1D6E4})$
parametrizes the upper sheet, so
$z_{\unicode[STIX]{x1D6E4}}=x_{\unicode[STIX]{x1D6E4}}+\text{i}y_{\unicode[STIX]{x1D6E4}}$
is a tangent. Hence by (2.5) the pseudo-velocity
$w^{\ast }-\unicode[STIX]{x1D707}z=q^{x}+\text{i}q^{y}$
is everywhere tangent. This represents the pseudo-velocity on the sheet; the limits on each side of the sheet differ by adding or subtracting
${\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{s}\boldsymbol{s}$
, where
$\boldsymbol{s}$
is a sheet unit tangent, subscript
$s$
indicates the derivative with respect to arc length; hence the one-sided limits of pseudo-velocity are also tangential. The condition
$\boldsymbol{q}\boldsymbol{\cdot }\boldsymbol{n}=0$
, i.e.
$\boldsymbol{v}\boldsymbol{\cdot }\boldsymbol{n}=\unicode[STIX]{x1D707}\boldsymbol{x}\boldsymbol{\cdot }\boldsymbol{n}$
, corresponds to
$\tilde{\boldsymbol{v}}\boldsymbol{\cdot }\boldsymbol{n}=\unicode[STIX]{x1D70E}$
for general unsteady vortex sheets,
$\unicode[STIX]{x1D70E}$
being the normal speed of the sheet.
When considering flows in the entire two-dimensional plane we assume the velocity
$w$
is bounded or at least
$o(|z|)$
near infinity (as
$|z|\rightarrow \infty$
), so that velocity is uniquely determined by vorticity. If so, then the term
$-\unicode[STIX]{x1D707}z$
in (2.5) dominates
$w$
near infinity, leading to an approximation
$$\begin{eqnarray}\displaystyle 0=(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}z_{\unicode[STIX]{x1D6E4}}-\unicode[STIX]{x1D707}z, & & \displaystyle\end{eqnarray}$$
which has a simple solution
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6E4}=\text{const.}\cdot |z|^{2-1/\unicode[STIX]{x1D707}},\quad \frac{z}{|z|}=\text{const}. & & \displaystyle\end{eqnarray}$$
Hence, near infinity, the vortex sheets become asymptotic to half-lines which are exactly those in the initial data (figure 1
a); the parameter
$\unicode[STIX]{x1D719}_{\infty }$
is determined by the initial data as half the angle between the two sheets in each pair.
2.2 Horizontal velocity approximation
2.2.1 Cusp contributions
We consider cusp solutions of (2.5). A key idea (see figure 4
a) is that the velocity induced in a point
$z_{0}$
by a not too close point
$z^{\prime }$
is almost completely cancelled by the mirror image
${z^{\prime }}^{\ast }$
. More precisely, an infinitesimal segment
$\text{d}\unicode[STIX]{x1D6E4}^{\prime }$
at
$z^{\prime }$
on the upper sheet induces
$$\begin{eqnarray}\displaystyle \frac{\text{d}\unicode[STIX]{x1D6E4}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}(z_{0}-x^{\prime }-\text{i}y^{\prime })}, & & \displaystyle\end{eqnarray}$$
whereas its mirror image induces
$$\begin{eqnarray}\displaystyle \frac{-\text{d}\unicode[STIX]{x1D6E4}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}(z_{0}-x^{\prime }+\text{i}y^{\prime })}; & & \displaystyle\end{eqnarray}$$
the sum is
$$\begin{eqnarray}\displaystyle \text{d}w=\frac{\text{d}\unicode[STIX]{x1D6E4}^{\prime }}{\unicode[STIX]{x03C0}}\frac{y^{\prime }}{(z_{0}-x^{\prime })^{2}+{y^{\prime }}^{2}}. & & \displaystyle\end{eqnarray}$$
The value of
$y^{\prime }$
is very small near the cusp, so the sum is much smaller than the summands unless
$z_{0}$
is very close to
$z^{\prime }$
.
We emphasize the importance of this cancellation: the velocity integral, an operator that is generally non-local and thus complicated, is in the cusp case dominated by its local part. Of course ‘local’ nurtures hopes of ‘differential’, and indeed this convenient reduction will soon become apparent.
In addition, for
$z_{0}$
on the upper sheet there is cancellation between the ‘near’ parts to the left and to the right as well (see figure 4
b): for positive
$\unicode[STIX]{x1D6E4}_{x}$
the left part induces almost upward velocities and the right one almost downward ones. So the dominant contribution to velocity is from the nearby mirror image.
Figure 5. (a) Sheet tangents become almost horizontal near the cusp; (b) velocity field for exactly horizontal sheets with constant velocity jumps
$\pm \unicode[STIX]{x1D6E4}_{x}$
.
To calculate its effect the following argument is convenient: in near-cusp points
$x+\text{i}y$
the sheet tangent is almost horizontal (figure 5
a); hence we approximate the sheet and its mirror image by two straight lines with constant velocity jump
$\unicode[STIX]{x1D6E4}_{x}(x)$
(figure 5
b). Such an approximation has velocity
$\overline{\boldsymbol{v}}=(0,0)$
above the upper sheet (and below the lower one),
$\text{}\underline{\boldsymbol{v}}=(\unicode[STIX]{x1D6E4}_{x},0)$
between the two sheets; on the upper sheet the principal value is the arithmetic average
$\boldsymbol{v}=({\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{x},0)$
. For the real (horizontal) part of (2.5),
$$\begin{eqnarray}\displaystyle 0 & = & \displaystyle v^{x}-\unicode[STIX]{x1D707}x+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}x_{\unicode[STIX]{x1D6E4}}\nonumber\\ \displaystyle & {\approx} & \displaystyle {\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{x}-\unicode[STIX]{x1D707}x+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}x_{\unicode[STIX]{x1D6E4}}.\end{eqnarray}$$
In one case our approximations are appropriate not only near the cusp point but everywhere: let
$2\unicode[STIX]{x1D719}_{\infty }$
be the angle the straight-line asymptotes of vortex sheet and mirror image enclose near infinity (figure 1
a). If the infinity angle
$\unicode[STIX]{x1D719}_{\infty }$
is small, then we may expect that the sheet tangents are nearly horizontal everywhere.
2.2.2 Non-cusp contributions; symmetry
Of course we have neglected that there are velocity contributions from the outer (non-cusp) part of the sheets, as well as those of other sheets, vortex patches or any other form of vorticity present in the two-dimensional similarity plane. Since they are uniformly distant from the cusp we may model the velocity induced by them as a complex-analytic function
$e=e(z)$
. We use its Taylor expansion
$$\begin{eqnarray}\displaystyle e(z)=e_{0}+e_{1}(z-x_{\ast })+o(|z-x_{\ast }|), & & \displaystyle\end{eqnarray}$$
where
$x_{\ast }$
is the cusp location which is real (horizontal axis); the coefficients
$e_{0},e_{1}$
are likewise real since velocity is horizontal at the wall by the slip condition. Expansions beyond
$e_{1}$
cannot always be justified (Elling Reference Elling2019b
, § 4.1). The constant
$e_{0}$
does not have a direct influence on the cusp shape; it merely corresponds to a constant shift in the similarity plane. The strain coefficient
$e_{1}$
corresponds to a saddle flow
$\boldsymbol{v}=(e_{1}(x-x_{\ast }),-e_{1}y)$
; for
$e_{1}>0$
this flow is expanding along the symmetry axis of the cusp (horizontal axis), compressing in the perpendicular direction. Equation (2.11) is replaced by the more accurate model
$$\begin{eqnarray}\displaystyle 0={\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{x}+\underbrace{e_{0}+e_{1}\unicode[STIX]{x1D709}}_{{\approx}\operatorname{Re}e}-\unicode[STIX]{x1D707}x+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}x_{\unicode[STIX]{x1D6E4}}. & & \displaystyle\end{eqnarray}$$
In a few special cases we may make additional assumptions. For example, for small infinity angle
$\unicode[STIX]{x1D719}_{\infty }$
, and without non-cusp flow features present, we may assume that
$e$
becomes small (possibly in contrast to
${\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{x}$
). Another special case is
$N$
-fold symmetry (see figure 1): instead of a single pair of equal-strength opposite-circulation sheets consider
$N$
pairs, with equal angles
$2\unicode[STIX]{x03C0}/N$
between the symmetry axes of adjacent pairs and cusps in the origin (as suggested by numerical experiments). Then the BR equation is
$$\begin{eqnarray}\displaystyle & \displaystyle 0=w^{\ast }-\unicode[STIX]{x1D707}z+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}z_{\unicode[STIX]{x1D6E4}}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle w=\text{p.v.}\int \mathop{\sum }_{n=0}^{N-1}\left(\frac{1}{z(\unicode[STIX]{x1D6E4})-\text{e}^{2\unicode[STIX]{x03C0}\text{i}n/N}z(\unicode[STIX]{x1D6E4}^{\prime })}-\frac{1}{z(\unicode[STIX]{x1D6E4})-(\text{e}^{2\unicode[STIX]{x03C0}\text{i}n/N}z(\unicode[STIX]{x1D6E4}^{\prime }))^{\ast }}\right)\frac{\text{d}\unicode[STIX]{x1D6E4}^{\prime }}{2\unicode[STIX]{x03C0}\text{i}}, & \displaystyle\end{eqnarray}$$
where
$n$
indexes the pair and
$z=z(\unicode[STIX]{x1D6E4})$
parametrizes the upper sheet of pair
$n=0$
; it is sufficient to solve the equation on one sheet, then the corresponding equation for the image and every other pair holds by symmetry. The effect on each sheet pair from near-cusp parts of other pairs can again be neglected due to the aforementioned strong
$\pm$
cancellation. The Taylor expansion of an
$N$
-symmetric
$e$
vanishes up to and including
$z^{N-1}$
. For
$N=2$
that means
$e_{0}=0$
; for
$N\geqslant 3$
we additionally have
$e_{1}=0$
. Symmetry is a convenient device for suppressing
$e_{0},e_{1}$
so that other influences shaping the cusp, such as dependence on
$\unicode[STIX]{x1D707}$
or
$\unicode[STIX]{x1D719}_{\infty }$
, can be studied in isolation.
Symmetry also helps in avoiding problems at infinity: the integral in the self-similar Birkhoff–Rott equation (2.5) behaves like
$\int |z|^{-1/\unicode[STIX]{x1D707}}\,\text{d}|z|$
near infinity, which generally diverges for any
$\unicode[STIX]{x1D707}\geqslant 1$
. In special cases it does converge; for instance with
$N\geqslant 2$
symmetry the antipodal vortex sheets cancel out sufficiently to improve the integrand to
${\sim}\!|z|^{-1-1/\unicode[STIX]{x1D707}}$
, convergent for all
$\unicode[STIX]{x1D707}<\infty$
. Even so, infinity becomes dominant like
$O(\unicode[STIX]{x1D707})$
as
$\unicode[STIX]{x1D707}{\nearrow}\infty$
, but for
$N\geqslant 3$
a further cancellation avoids even that effect.
2.2.3 Conservation arguments
Equation (2.13) is not only a differential equation but an ordinary one; we have modelled
$v^{x}$
as a function of
$\unicode[STIX]{x1D709}$
alone. The ODE is not difficult to solve, but we minimize effort on irrelevant solutions by first considering conservation of mass;
$\unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{v}=0$
has a self-similar form
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{q}=\unicode[STIX]{x1D735}\boldsymbol{\cdot }(\boldsymbol{v}-\unicode[STIX]{x1D707}\boldsymbol{x})=-2\unicode[STIX]{x1D707}, & & \displaystyle\end{eqnarray}$$
where
$2$
comes from the number of dimensions. Integrate this over a small ‘cusp triangle’ (formed in figure 6 by solid and dashed vortex sheet and dashed-dotted ‘mouth’) from the cusp
$x_{\ast }$
to some
$x>x_{\ast }$
: since the normal pseudo-velocity
$\boldsymbol{q}\boldsymbol{\cdot }\boldsymbol{n}$
is zero at each side of a vortex sheet, only the ‘mouth’ boundary term remains,
$$\begin{eqnarray}\displaystyle \int _{-y(\unicode[STIX]{x1D709})}^{y(\unicode[STIX]{x1D709})}q^{x}\,\text{d}s=-2\unicode[STIX]{x1D707}\cdot \text{cusp triangle area}. & & \displaystyle\end{eqnarray}$$
On the left we use
$q^{x}\approx \text{}\underline{q}^{x}(\unicode[STIX]{x1D709})$
. The right-hand side is always negative, but for a ‘cusp’ the area is reasonably estimated as less than the area of the straightened-side triangle (dotted in figure 6) which is
$\unicode[STIX]{x1D709}y(\unicode[STIX]{x1D709})$
. Hence
$$\begin{eqnarray}\displaystyle 0>\text{}\underline{q}^{x}(\unicode[STIX]{x1D709})\geqslant -\unicode[STIX]{x1D707}\unicode[STIX]{x1D709}. & & \displaystyle\end{eqnarray}$$
Velocity
$v^{x}$
must be a constant plus
$O(\unicode[STIX]{x1D709})$
, and the same for
$\unicode[STIX]{x1D6E4}_{x}$
, which differs from
$v^{x}$
by a smooth function
$\operatorname{Re}e$
.
Figure 6. Mass conservation implies negative
$q^{x}$
averages over the right side of cusp triangles.
2.2.4 Constraints
Multiply both sides of (2.13) by
$\unicode[STIX]{x1D6E4}_{x}$
to obtain
$$\begin{eqnarray}\displaystyle 0={\textstyle \frac{1}{2}}\unicode[STIX]{x1D6E4}_{x}^{2}+(e_{0}+e_{1}\unicode[STIX]{x1D709}-\unicode[STIX]{x1D707}x)\unicode[STIX]{x1D6E4}_{x}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}; & & \displaystyle\end{eqnarray}$$
completing the square for
$\unicode[STIX]{x1D6E4}_{x}$
yields
$$\begin{eqnarray}\displaystyle 0={\textstyle \frac{1}{2}}(\unicode[STIX]{x1D6E4}_{x}+e_{0}+e_{1}\unicode[STIX]{x1D709}-\unicode[STIX]{x1D707}x)^{2}-{\textstyle \frac{1}{2}}(e_{0}+e_{1}\unicode[STIX]{x1D709}-\unicode[STIX]{x1D707}x)^{2}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}. & & \displaystyle\end{eqnarray}$$
The first parenthesis is recognized to be
$\text{}\underline{q}^{x}$
since (recall figure 5
b) the self-induced velocity is
$\unicode[STIX]{x1D6E4}_{x}$
between the sheets, so
$$\begin{eqnarray}\displaystyle \text{}\underline{q}^{x}=\text{}\underline{v^{x}}-\unicode[STIX]{x1D707}x=\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}+e_{0}+e_{1}\unicode[STIX]{x1D709}-\unicode[STIX]{x1D707}\underbrace{x}_{=x_{\ast }+\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$
So we obtain
$$\begin{eqnarray}\displaystyle 0={\textstyle \frac{1}{2}}(\text{}\underline{q}^{x})^{2}-{\textstyle \frac{1}{2}}(e_{0}+e_{1}\unicode[STIX]{x1D709}-\unicode[STIX]{x1D707}x)^{2}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}. & & \displaystyle\end{eqnarray}$$
Mass conservation was shown to require
$\text{}\underline{q}^{x}{\nearrow}0$
as we approach the cusp (
$\unicode[STIX]{x1D709}{\searrow}0$
), which implies (for
$\ast$
indicating values in the cusp limit)
$$\begin{eqnarray}\displaystyle (2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}_{\ast }={\textstyle \frac{1}{2}}(e_{0}-\unicode[STIX]{x1D707}x_{\ast })^{2}\geqslant 0. & & \displaystyle\end{eqnarray}$$
Hence
$\unicode[STIX]{x1D6E4}_{\ast }=0$
or
$\unicode[STIX]{x1D6E4}_{\ast }>0$
; the latter appears to happen only in borderline cases which are discussed in Elling (Reference Elling2019b
, § 4.2). Consider
$\unicode[STIX]{x1D6E4}_{\ast }=0$
from now on.
2.2.5 The
$x$
ODE solutions in case
$\unicode[STIX]{x1D6E4}_{\ast }=0$
By (2.23)
$\unicode[STIX]{x1D6E4}_{\ast }=0$
requires
$$\begin{eqnarray}\displaystyle e_{0}=\unicode[STIX]{x1D707}x_{\ast }, & & \displaystyle\end{eqnarray}$$
so (2.21) and (2.22) simplify to
$$\begin{eqnarray}\displaystyle \text{}\underline{q}^{x}=\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}+(e_{1}-\unicode[STIX]{x1D707})\unicode[STIX]{x1D709}\quad [\quad \Leftrightarrow \quad \unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}=\text{}\underline{q}^{x}+(\unicode[STIX]{x1D707}-e_{1})\unicode[STIX]{x1D709}\quad ] & & \displaystyle\end{eqnarray}$$
and
$$\begin{eqnarray}\displaystyle 0={\textstyle \frac{1}{2}}(\text{}\underline{q}^{x})^{2}-{\textstyle \frac{1}{2}}((e_{1}-\unicode[STIX]{x1D707})\unicode[STIX]{x1D709})^{2}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}. & & \displaystyle\end{eqnarray}$$
After
$\unicode[STIX]{x2202}_{\unicode[STIX]{x1D709}}$
of the latter, producing
$$\begin{eqnarray}\displaystyle 0=\text{}\underline{q}^{x}\text{}\underline{q}_{\unicode[STIX]{x1D709}}^{x}-(\unicode[STIX]{x1D707}-e_{1})^{2}\unicode[STIX]{x1D709}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}, & & \displaystyle\end{eqnarray}$$
the former substitutes the final
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}$
to yield
$$\begin{eqnarray}\displaystyle 0=\text{}\underline{q}^{x}\text{}\underline{q}_{\unicode[STIX]{x1D709}}^{x}+(\unicode[STIX]{x1D707}-1+e_{1})(\unicode[STIX]{x1D707}-e_{1})\unicode[STIX]{x1D709}+(2\unicode[STIX]{x1D707}-1)\text{}\underline{q}^{x}. & & \displaystyle\end{eqnarray}$$
This ODE turns autonomous if we eliminate the
$\unicode[STIX]{x1D709}$
dependence by the substitution
$$\begin{eqnarray}\displaystyle \text{}\underline{q}^{x}(\unicode[STIX]{x1D709})=\unicode[STIX]{x1D709}Q(\unicode[STIX]{x1D709}),\quad \unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}(\unicode[STIX]{x1D709})=\unicode[STIX]{x1D709}G(\unicode[STIX]{x1D709}), & & \displaystyle\end{eqnarray}$$
(which is also natural in light of
$\text{}\underline{q}^{x}=O(\unicode[STIX]{x1D709})$
) to obtain
$$\begin{eqnarray}\displaystyle 0=\unicode[STIX]{x1D709}^{2}QQ_{\unicode[STIX]{x1D709}}+\unicode[STIX]{x1D709}Q^{2}+(\unicode[STIX]{x1D707}-1+e_{1})(\unicode[STIX]{x1D707}-e_{1})\unicode[STIX]{x1D709}+(2\unicode[STIX]{x1D707}-1)\unicode[STIX]{x1D709}Q. & & \displaystyle\end{eqnarray}$$
Divide by
$\unicode[STIX]{x1D709}>0$
and factor using
$(\unicode[STIX]{x1D707}-1+e_{1})+(\unicode[STIX]{x1D707}-e_{1})=2\unicode[STIX]{x1D707}-1$
,
$$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}Q}{\unicode[STIX]{x2202}(-\log \unicode[STIX]{x1D709})}=\frac{-(Q+\unicode[STIX]{x1D707}-e_{1})(Q+\unicode[STIX]{x1D707}+e_{1}-1)}{-Q}. & & \displaystyle\end{eqnarray}$$
Our conservation-of-mass arguments showed that
$Q=\text{}\underline{q}^{x}/\unicode[STIX]{x1D709}$
is negative but bounded in the cusp limit
$-\log \unicode[STIX]{x1D709}\rightarrow +\infty$
(see (2.18)). This is possible only if
$Q$
converges to a non-positive root of the numerator. The two roots are
$$\begin{eqnarray}\displaystyle & \displaystyle Q=1-e_{1}-\unicode[STIX]{x1D707}\quad \overset{\text{(2.25)}}{{\rightsquigarrow}}\quad G=1-2e_{1}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle Q=e_{1}-\unicode[STIX]{x1D707}\quad \overset{\text{(2.25)}}{{\rightsquigarrow}}\quad G=0. & \displaystyle\end{eqnarray}$$
The two roots are equal if and only if
$e_{1}={\textstyle \frac{1}{2}}$
, a borderline case we ignore to keep the discussion concise and avoid
$\log \unicode[STIX]{x1D709}$
terms.
First, consider solutions that converge to one of the roots from one side. Then, that root must be asymptotically stable on that side. The right-hand side of (2.31) has positive denominator for
$Q<0$
whereas the leading
$-Q^{2}$
part of the quadratic numerator makes it positive between the roots, negative outside. Hence if the roots are distinct, the smaller root is either unstable on both sides or not negative and hence irrelevant, so we may ignore it. On the other hand, if the larger root is negative, then it is asymptotically stable on both sides. The corresponding solutions will indeed turn out to lead to meaningful vortex cusps.
If
$e_{1}<{\textstyle \frac{1}{2}}$
, then the
$Q=1-e_{1}-\unicode[STIX]{x1D707}$
root is the larger root. It is negative if and only if
$e_{1}>1-\unicode[STIX]{x1D707}$
. We only consider
$\unicode[STIX]{x1D707}>{\textstyle \frac{1}{2}}$
(see the introduction), so
$1-\unicode[STIX]{x1D707}<{\textstyle \frac{1}{2}}$
, meaning the range
$]1-\unicode[STIX]{x1D707},{\textstyle \frac{1}{2}} [$
for
$e_{1}$
is never empty. However, the constraint
$1-\unicode[STIX]{x1D707}<e_{1}$
is significant, especially for smaller
$\unicode[STIX]{x1D707}$
; for
$\unicode[STIX]{x1D707}<1$
it does not permit cusps with
$e_{1}$
being arbitrarily small. The root corresponds to
$G=1-2e_{1}+\cdots \,$
, so we obtain solutions
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}=(1-2e_{1})\unicode[STIX]{x1D709}+o(\unicode[STIX]{x1D709}),\quad \unicode[STIX]{x1D6E4}=({\textstyle \frac{1}{2}}-e_{1})\unicode[STIX]{x1D709}^{2}+o(\unicode[STIX]{x1D709}^{2}). & & \displaystyle\end{eqnarray}$$
If on the other hand
$e_{1}>{\textstyle \frac{1}{2}}$
, then
$Q=e_{1}-\unicode[STIX]{x1D707}$
is the larger root. It is negative if and only if
$e_{1}<\unicode[STIX]{x1D707}$
; again the range
$]{\textstyle \frac{1}{2}},\unicode[STIX]{x1D707}[$
is non-empty, but especially for smaller
$\unicode[STIX]{x1D707}$
, a significant constraint. The root corresponds to
$G=0+\cdots \,$
so the leading term of
$\unicode[STIX]{x1D6E4}$
has to be determined by linearizing (2.31) there: passing to
$$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}G}{\unicode[STIX]{x2202}(-\log \unicode[STIX]{x1D709})}\overset{\text{(2.31)}}{=}\frac{-G(G+2e_{1}-1)}{-G+\unicode[STIX]{x1D707}-e_{1}} & & \displaystyle\end{eqnarray}$$
the linearization at
$G=0$
is
$$\begin{eqnarray}\displaystyle \frac{\unicode[STIX]{x2202}G}{\unicode[STIX]{x2202}(-\log \unicode[STIX]{x1D709})}=-\frac{2e_{1}-1}{\unicode[STIX]{x1D707}-e_{1}}G, & & \displaystyle\end{eqnarray}$$
with solutions
$G\sim \unicode[STIX]{x1D709}^{(2e_{1}-1)/(\unicode[STIX]{x1D707}-e_{1})}$
.
In summary,
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}=C\unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FD}-1}+\cdots ,\quad \unicode[STIX]{x1D6E4}=\frac{C}{\unicode[STIX]{x1D6FD}}\unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FD}}+\cdots \,, & & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle C=\left\{\begin{array}{@{}ll@{}}1-2e_{1},\quad & e_{1}<{\textstyle \frac{1}{2}}\\ \text{free},\quad & e_{1}>{\textstyle \frac{1}{2}}.\end{array}\right. & & \displaystyle\end{eqnarray}$$
Finally, consider solutions that equal one of the roots, so that asymptotic stability is not necessary. The case
$Q=e_{1}-\unicode[STIX]{x1D707}$
corresponds to
$\unicode[STIX]{x1D6E4}_{x}=0$
, i.e. the trivial case of no sheets at all, which we ignore. The case
$Q=1-e_{1}-\unicode[STIX]{x1D707}$
for
$e_{1}<{\textstyle \frac{1}{2}}$
is a special case of the discussion above; for
$e_{1}>{\textstyle \frac{1}{2}}$
it does not produce cusps, as will become apparent at the end of the following section.
2.3 Vertical velocity modelling
Given the relationship between
$\unicode[STIX]{x1D709}$
and
$\unicode[STIX]{x1D6E4}$
, it remains to derive a model for the upper-sheet height
$y=y(\unicode[STIX]{x1D709})$
. Two natural and simple but flawed approximations are discussed elsewhere (Elling Reference Elling2019b
, § 2.3); here, we focus on the correct approach: the requirement of mass conservation only leaves one possible choice. Consider again the cusp triangle argument (figure 6),
$$\begin{eqnarray}\displaystyle -2\unicode[STIX]{x1D707}\cdot \text{area}=2y\text{}\underline{q}^{x}. & & \displaystyle\end{eqnarray}$$
Take
$\unicode[STIX]{x2202}_{\unicode[STIX]{x1D709}}$
, corresponding to moving the cusp triangle ‘mouth’,
$$\begin{eqnarray}\displaystyle -2\unicode[STIX]{x1D707}\cdot 2y=2y_{\unicode[STIX]{x1D709}}\text{}\underline{q}^{x}+2y\text{}\underline{q}_{\unicode[STIX]{x1D709}}^{x}. & & \displaystyle\end{eqnarray}$$
Solve for
$$\begin{eqnarray}\displaystyle \frac{y_{\unicode[STIX]{x1D709}}}{y}=\frac{2\unicode[STIX]{x1D707}+\text{}\underline{q}_{\unicode[STIX]{x1D709}}^{x}}{-\text{}\underline{q}^{x}}. & & \displaystyle\end{eqnarray}$$
Substitute (2.25),
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D709}\frac{y_{\unicode[STIX]{x1D709}}}{y}=\frac{\unicode[STIX]{x1D707}+e_{1}+\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}\unicode[STIX]{x1D709}}}{\unicode[STIX]{x1D707}-e_{1}-\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}/\unicode[STIX]{x1D709}}. & & \displaystyle\end{eqnarray}$$
Equation (2.37) shows
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D709}\frac{y_{\unicode[STIX]{x1D709}}}{y}=\unicode[STIX]{x1D6FC}+o(1), & & \displaystyle\end{eqnarray}$$
with solution
$y=\text{const.}\cdot \unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FC}}$
, where
We emphasize that
$e_{1}$
is not a ‘fudge factor’ that can be used to realize any observed cusp exponent
$\unicode[STIX]{x1D6FC}$
; e.g. in case of
$N\geqslant 3$
symmetry
$e_{1}$
is zero so that the prediction (2.45) depends only on
$\unicode[STIX]{x1D707}$
. In non-symmetric numerics
$e_{1}$
can be measured independently; the formula (2.38) for
$\unicode[STIX]{x1D6FD}$
(also measurable) is likewise dependent on
$e_{1}$
and
$\unicode[STIX]{x1D707}$
, providing a second constraint for the two. Here,
$e_{1}$
cannot be modelled locally since it depends on the vorticity in every part of the similarity plane.
At this point we may collect the loose end left in § 2.2.5 where we ignored the constant solution
$Q=1-e_{1}-\unicode[STIX]{x1D707}$
, i.e.
$\unicode[STIX]{x1D6E4}_{x}=(1-2e_{1})x$
, in case
$e_{1}>{\textstyle \frac{1}{2}}$
. Repeating our conservation of mass argument we obtain
$$\begin{eqnarray}\displaystyle \unicode[STIX]{x1D6FC}=\frac{\unicode[STIX]{x1D707}+1-e_{1}}{\unicode[STIX]{x1D707}+e_{1}-1}\overset{e_{1}>{\textstyle \frac{1}{2}}}{{\leqslant}}1. & & \displaystyle\end{eqnarray}$$
But a proper
$y\sim x^{\unicode[STIX]{x1D6FC}}$
cusp, in particular our small-slope assumption
$|y_{x}|\ll 1$
, requires
$\unicode[STIX]{x1D6FC}>1$
. Although this
$Q$
represents an appropriate solution of the
$x$
ODE, it leads to an unusable solution of the
$y$
ODE. Hence we were justified in rejecting it.
3 Numerical validation
3.1 Methods
To test our model we compare its predictions to numerical data obtained from discretizing the Birkhoff–Rott equation (2.15) itself. Besides general methods for unsteady vortex sheet rollup or other free-surface evolutions (e.g. Meiron, Baker & Orszag Reference Meiron, Baker and Orszag1982; Krasny Reference Krasny1991; Nitsche Reference Nitsche2001), numerics adapted to self-similar flow using point vortices were successfully used to compute algebraic vortex spirals in Pullin (Reference Pullin1978, Reference Pullin and Caflisch1989). The numerical scheme parametrizes the upper sheet by a function
$z$
of
$\unicode[STIX]{x1D6E4}$
, solving for variable point vortices
$z_{i}$
approximating
$z(\unicode[STIX]{x1D6E4}_{i})$
in finitely many
$\unicode[STIX]{x1D6E4}_{i}$
. Evaluation points for the BR equation are midway between adjacent
$z_{i}$
to mitigate the strong singularity of the point vortex velocity field.
Figure 7. Point vortex numerics: velocities are singular near each
$\pm$
pair, near-zero elsewhere.
In our experience this approach is quite effective for vortex spirals, but for cusps or similar cases it is impractical due to the cancellation between positive and negative vorticity (figure 7): unless a prohibitively expensive number of vortex pairs is used, near the cusp the horizontal distance between vortex pairs must be much larger than the vertical distance between vortex and mirror image in each pair. Then the velocity field would be highly singular in the vicinity of a pair, but nearly zero elsewhere due to cancellation – clearly no resemblance to the correct field for smooth sheets. Inevitably approximation by higher-degree polynomials must be used instead.
The sheet segment from
$z_{i}$
to
$z_{i+1}$
is approximated by a polynomial interpolating these points and a few adjacent ones. The precise choice of
$z_{i}$
is crucial for point vortex methods, but not for the higher-order methods used here. The BR integral
$\int (z-z^{\prime })^{-1}\,\text{d}\unicode[STIX]{x1D6E4}^{\prime }$
can be done by partial fractions, requiring complex logarithms and some root finding, e.g. using quadratic or Cardano formulas. (In this article we are focused on modelling correctness rather than numerical speed.) Equation (2.15) is discretized by evaluation in the arithmetic averages
$\unicode[STIX]{x1D6E4}_{i+1/2}$
of
$\unicode[STIX]{x1D6E4}_{i},\unicode[STIX]{x1D6E4}_{i+1}$
, which produces finitely many nonlinear algebraic equations with complex values depending on complex variables
$z_{i}$
. Newton iteration is applied to the system, using direct methods for linear solves with the exact derivative matrix. Cubic cost growth limits the number of vortices to a few thousand, which is sufficient for our purposes. Spot checks with fine horizontal resolution showed that no significant high-frequency features were suppressed.
Figure 8. Case
$\unicode[STIX]{x1D707}=1.3$
with
$N=2$
symmetry. (a)
$\unicode[STIX]{x1D719}_{\infty }=10^{\circ }$
a cusp forms (arrows
$\boldsymbol{v}$
). (b)
$\unicode[STIX]{x1D719}_{\infty }=25.5^{\circ }$
, near the limit angle (arrows
$\boldsymbol{q}$
; nearly zero between the sheets where they meet). (c)
$\unicode[STIX]{x1D719}_{\infty }=35^{\circ }$
; instead of a cusp a tiny jet has formed (arrows
$\boldsymbol{v}$
); (d) jet detail (arrows
$\boldsymbol{q}$
; saddle points at jet entrance and exit). (e)
$\unicode[STIX]{x1D719}_{\infty }=48^{\circ }$
jet, still narrow (arrows
$\boldsymbol{v}$
). (f)
$\unicode[STIX]{x1D719}_{\infty }=77^{\circ }$
jet (arrows
$\boldsymbol{v}$
).
The exact sheet
$\unicode[STIX]{x1D6E4}$
ranges from
$0$
(cusp) to
$\infty$
; the approximation uses inner and outer cutoffs. At the outer cutoff several fixed
$z$
are calculated by the straight-sheet asymptotes (2.7). The outer cutoff must be taken large enough to avoid large errors, especially for
$N\leqslant 2$
and large
$\unicode[STIX]{x1D707}$
. At the inner cutoff, for the innermost segment
$z_{0},z_{1}$
the polynomial degree was limited to quadratic using
$z_{2}$
as third knot.
Like Pullin’s point vortex calculations ours are very sensitive due to the inherent instability of vortex sheets. Sheet and image are prone to intersecting each other or to separating at the inter cutoff during iteration, causing uncontrollable oscillations. It is important to start with a good initial approximation and to limit the iteration step size until a good basin of attraction is reached. Initial approximations can be obtained by solving our ODE models (2.13) and (2.43). The ODE initial conditions are imposed at some large
$\unicode[STIX]{x1D6E4}$
using the asymptotics (2.7) that are valid near infinity; the first constant is the sheet strength, which can be chosen arbitrarily, while the second constant
$z/|z|$
corresponds to
$\unicode[STIX]{x1D719}_{\infty }$
, half the angle between the two sheets in each pair near infinity (see figure 1). Recall that infinity in the similarity plane corresponds to the initial data in
$(t,x)$
coordinates. The ODE model solutions are especially accurate for small
$\unicode[STIX]{x1D719}_{\infty }$
; once the numerics have converged
$\unicode[STIX]{x1D719}_{\infty }$
can gradually be increased to larger values. For additional details, such as calculation of spiral-ended sheets, see (Elling Reference Elling2019b
, § 3).
3.2 Results
Vector fields in all diagrams use a truncated logarithmic scaling to avoid very large or very small arrows. Numerical spirals shown in this paper may use an inner cutoff to a central point vortex, leaving an empty core in the figures, as in figure 8(d) physical spirals without cutoff look like figure 2(c), with a filled-in core.
Input parameters for each problem are similarity exponent
$\unicode[STIX]{x1D707}$
, number of sheet pairs
$N$
and the angle
$\unicode[STIX]{x1D719}_{\infty }$
between sheets in each pair near infinity. We do not specify the sheet strength since it can be scaled arbitrarily by scaling time
$t$
and space
$x$
. Vortex spacing, cutoffs and degree are chosen to yield results that do not change significantly upon further refinement, unless otherwise noted.
Our numerical results (see figure 8
a) show that vortex cusps exist under some conditions, namely if
$\unicode[STIX]{x1D707}$
is sufficiently large, and if the data are not too ‘large’, which means in the special case of global cusps that
$\unicode[STIX]{x1D719}_{\infty }$
is sufficiently small, below a limit angle that depends on
$\unicode[STIX]{x1D707}$
and other parameters. In cases where cusps do not appear spiral-ended jets can be observed instead.
3.2.1 Linear and higher degree
Figure 9. The
$N=3$
global cusp exponents
$\unicode[STIX]{x1D6FC}$
(vertical axis) from piecewise linear or quadratic numerics.
Cusp exponents were estimated by fixing two sufficiently distant points well inside the cusp region of the sheet pair, but not too close to the inner cutoff, to avoid spurious influences. Since for
$N=3$
we know the cusp location is the origin, fitting
$y=C_{1}x^{\unicode[STIX]{x1D6FC}}$
(or
$\unicode[STIX]{x1D6E4}=C_{2}x^{\unicode[STIX]{x1D6FD}}$
) to the data only requires taking logarithms and solving a simple linear equation for
$\log C_{i},\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D6FD}$
.
The lowest-degree linear interpolation of
$z_{i},z_{i+1}$
, consistently produces cusps with incorrect exponent (dotted curve in figure 9), a value that corresponds to the exponent from an oversimplified conservation-violating model (Elling (Reference Elling2019b
, § 2.3.2), dotted curve in figure 9), as confirmed by numerical mass balance checks. Convergence/stability especially for small
$\unicode[STIX]{x1D707}$
is poor, contrary to the default expectation that lower-order numerics are more robust. Linear approximation is unable to resolve the term
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}\unicode[STIX]{x1D709}}$
in (2.43).
Quadratic reconstructions, as well as higher-degree ones, are more robust and yield consistent results, in particular the exponent (2.45) (solid curve in figure 9, closely matching the ‘
$\times$
’ for quadratic numerics). The match is good except near
$\unicode[STIX]{x1D707}=1$
where large cusp exponent and cusp–jet transition start to take effect.
The logarithmic singularities at corners
$z_{i}$
did not appear to hurt accuracy, similar to (Pullin Reference Pullin1978, Reference Pullin and Caflisch1989) where the far more severe
$1/r$
singularities of point vortices still permitted accurate computation of vortex spirals. Enforcing continuous tangents or higher derivatives caused instability without improving accuracy; e.g. natural cubic splines are well known to oscillate near low-regularity features like the fork point in figure 8(b), causing sheet and image to intersect.
Quadratic or higher-degree approximation is only necessary for the segment containing the evaluation point and for its mirror image; linear or even points are sufficient elsewhere. This does not apply to more complex flows where a third segment may be close to the evaluation point. We use either quadratic or quadratic-nearby-linear-elsewhere reconstruction for the rest of this paper; experiments did not show significant differences between them or cubic or higher degrees.
3.2.2 Limit angle
The limit angle for
$N=3$
global cusps is shown in figure 10. Clearly, below
$\unicode[STIX]{x1D707}=1.3$
the limit angles become rather restrictive, converging to
$0$
as
$\unicode[STIX]{x1D707}{\searrow}1$
. For
$\unicode[STIX]{x1D707}\leqslant 1.05$
, near the left side, large cusp exponents cause unreliable calculations. The limit angle diagram for
$N=2$
and other cases is similar; in those cases
$e_{1}\rightarrow 0$
along with
$\unicode[STIX]{x1D707}{\searrow}1$
so that it is not possible to obtain cusps at or below
$\unicode[STIX]{x1D707}=1$
.
Figure 10. Upper limit on
$\unicode[STIX]{x1D719}_{\infty }$
for existence of
$N=3$
global cusps.
Figure 11. Case
$\unicode[STIX]{x1D707}=0.9$
,
$N=1$
(no symmetry), external
$w$
field adjusted to yield
$e_{1}\approx 0.3$
(arrows
$\boldsymbol{v}$
).
However, it is possible to see local cusps for
$\unicode[STIX]{x1D707}\leqslant 1$
if a large strain coefficient
$e_{1}$
is generated by other vorticity in the similarity plane or by externally imposed harmonic velocity fields. In figure 11 we add to the BR integral velocity an additional harmonic velocity field to obtain
$e_{1}\approx 0.3$
. It is possible to achieve cusps at
$\unicode[STIX]{x1D707}=0.9$
and lower, at fairly ‘large’ angles between the straight outer segments. The two straight segments on the right of figure 11 are fixed, only the curved parts solve the BR equations. (The pronounced corner between straight and curved segments was inconsequential.) Such local cusps can generally not be extended to infinity. Truncating the sheet at a low outer cutoff scale causes a complex-analytic change to the velocity near the cusp, so the effect can be modelled by adjusting our
$e$
. In addition, to stabilize the solution against outer-sheet changes, it is possible to subtract from all calculated
$w$
the value computed at a fixed point, say the origin. A constant added to
$w$
can in (2.15) be absorbed by a constant shift added to
$z$
, since every term except
$\unicode[STIX]{x1D707}z$
is either a difference of two
$z$
or a derivative of
$z$
. Hence only the location of the cusp is changed, but not its shape or circulation distribution.
3.2.3 The
$y$
,
$x$
exponents with strain
We also wish to demonstrate the effect of the ‘saddle coefficient’
$e_{1}$
on cusp exponents (figure 12). We use
$N=2$
, which permits
$e_{1}\neq 0$
while avoiding velocity integral divergence near infinity, and
$\unicode[STIX]{x1D719}_{\infty }=10^{\circ }$
, which is large enough to generate significant
$e_{1}$
but small enough to allow
$\unicode[STIX]{x1D707}$
near
$1$
. The coefficient
$e_{1}$
was estimated numerically: our approximation leading to (2.11) is that, in the velocity integral
$w$
in (2.15), contributions from near-cusp parts of sheets almost cancel so that
$w\approx e$
at points near the cusp but outside any sheet pair. Hence, the
$x$
derivative of the velocity integral in such points can be used to estimate the first-order Taylor coefficient
$e_{1}$
. The derivative can also be approximated by a finite difference in two distinct points that are not too close (high-precision arithmetic may be needed to avoid round-off errors).
Figure 12. Numerical cusp exponents
$\unicode[STIX]{x1D6FC}$
on the vertical axis in the presence of saddle strain.
We evaluate the cusp exponent
$\unicode[STIX]{x1D6FC}=(\unicode[STIX]{x1D707}+1-e_{1})/(\unicode[STIX]{x1D707}-1+e_{1})$
from (2.45) once with
$e_{1}=0$
(dotted curve in figure 12) and once with the numerically estimated
$e_{1}$
(solid curve); none of these
$e_{1}$
exceeded
${\textstyle \frac{1}{2}}$
. Numerically measured
$\unicode[STIX]{x1D6FC}$
fit the solid curve well, so we conclude that our modelling likely yielded the correct
$e_{1}$
dependence of the cusp exponent. The numerical cusp exponents fit the solid curve closely, except near the
$\unicode[STIX]{x1D707}{\searrow}1$
end, where large cusp exponents and transition to jetting cause now-familiar inaccuracy; the fit at large
$\unicode[STIX]{x1D707}$
can be improved further by taking larger outer cutoffs. This corresponds to
$e_{1}$
being more sensitive to outer parts of the sheets when
$\unicode[STIX]{x1D707}$
is large. To improve the fit from
$\unicode[STIX]{x1D707}=1.5$
to
$\unicode[STIX]{x1D707}=3$
it is sufficient to refine the discretization alone.
3.2.4 The
$\unicode[STIX]{x1D6E4}$
,
$x$
exponents with strain
Finally, in figure 13 we consider the exponent
$\unicode[STIX]{x1D6FD}$
in the
$\unicode[STIX]{x1D6E4}=C\unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FD}}+\cdots \,$
relationship in the cusp limit. Formula (2.38) distinguishes the cases
$e_{1}<{\textstyle \frac{1}{2}}$
and
$e_{1}>{\textstyle \frac{1}{2}}$
. The strain coefficient
$e_{1}$
induced by the non-cusp parts of the sheet increases monotonically as the infinity angle
$2\unicode[STIX]{x1D719}_{\infty }$
between the outer parts increases;
$e_{1}$
does exceed
${\textstyle \frac{1}{2}}$
if
$\unicode[STIX]{x1D707}$
is chosen sufficiently large (with
$\unicode[STIX]{x1D707}=2.5$
chosen for computing figure 13), allowing not only a larger angle but also using that, for
$\unicode[STIX]{x1D707}{\nearrow}\infty$
, the effects of the outer sheet in the velocity integral become ever stronger.
Figure 13. Solid curve: theoretical
$\unicode[STIX]{x1D6FD}$
(2.37); dotted curves: numerically estimated
$(e_{1},\unicode[STIX]{x1D6FD})$
.
The results, again using a numerically estimated
$e_{1}$
, show a very close match between formula (2.38) (solid curve) and the numerically determined
$\unicode[STIX]{x1D6FD}$
(dashed and dotted curves). Near
$e_{1}={\textstyle \frac{1}{2}}$
the results are almost
$0.02$
apart. Errors decrease rapidly at a distance from
$e_{1}={\textstyle \frac{1}{2}}$
as the inner cutoff is taken smaller and simultaneously the initial
$-\text{d}\unicode[STIX]{x1D709}/\unicode[STIX]{x1D709}$
step size is reduced (the dotted curve has inner cutoff improved by a factor
$1/10$
and resolution by
$1/2$
). But convergence is very slow near
$e_{1}\approx {\textstyle \frac{1}{2}}$
, which is natural: not only is it difficult to distinguish a power law exponent
$\unicode[STIX]{x1D6FD}=0.5$
from
$\unicode[STIX]{x1D6FD}=0.49$
numerically, but as
$e_{1}{\nearrow}{\textstyle \frac{1}{2}}$
the coefficient
$C$
in
$\unicode[STIX]{x1D6E4}=Cx^{0.5}+\cdots \,$
also vanishes (see (2.39)), meaning the leading term is dominated by the next one for all but astronomically small inner cutoffs.
We conclude that our formulas (2.38) and (2.45) for the
$\unicode[STIX]{x1D6E4},x$
and
$y,x$
exponents are probably the physically correct exponents, with close match between numerical data and heuristic derivations. We emphasize again that the
$y,x$
exponent is completely determined by the conservation of mass argument, so that the point of attack for criticism is confined to the
$\unicode[STIX]{x1D6E4},x$
exponent and the ODE model and simplifying assumptions that led to it.
3.2.5 Clockwise upper circulation
Figure 14. Spiral-ended jets for clockwise-circulation upper sheets (
$\unicode[STIX]{x1D707}=1.3$
,
$N=1$
,
$\unicode[STIX]{x1D719}_{\infty }=10^{\circ }$
, arrows
$\boldsymbol{q}$
).
Figure 15. Local cusp with clockwise upper right sheet circulation, made possible by strain coefficient
$e_{1}>{\textstyle \frac{1}{2}}$
; arrows
$\boldsymbol{v}$
.
For
$e_{1}<{\textstyle \frac{1}{2}}$
our
$x$
ODE discussion in § 2.2.5 yielded solutions
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}=({\textstyle \frac{1}{2}}-e_{1})\unicode[STIX]{x1D709}+\cdots \,$
; note
${\textstyle \frac{1}{2}}-e_{1}>0$
so that
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}>0$
near the cusp, corresponding to counterclockwise circulation on the upper sheet near cusps opening to the right. For clockwise circulation, jets rather than cusps will be observed, as in the local flow in figure 14 (
$\unicode[STIX]{x1D707}=1.3$
,
$N=1$
).
However, if
$e_{1}>{\textstyle \frac{1}{2}}$
, then § 2.2.5 permitted solutions
$\unicode[STIX]{x1D6E4}_{\unicode[STIX]{x1D709}}=C\unicode[STIX]{x1D709}^{\unicode[STIX]{x1D6FD}-1}+\cdots \,$
for
$\unicode[STIX]{x1D6FD}-1>1$
(see (2.38)), with arbitrary
$C$
. In particular, we may choose
$C<0$
, allowing clockwise circulation on the upper sheet at the cusp. Numerically it can be confirmed that such solutions are possible: figure 15 shows a local
$N=2$
symmetric cusp with
$\unicode[STIX]{x1D707}=1.5$
and
$e_{1}$
adjusted to
$0.7$
by adding an external harmonic velocity field
$w=cz$
with suitable constant
$c$
. Again the straight outer segments are fixed chosen data, whereas the curved inner parts are numerical Birkhoff–Rott solutions forming a cusp at the origin.
Generally, the case of global flows with clockwise upper sheet circulation and infinity angle
$\unicode[STIX]{x1D719}_{\infty }$
is equivalent to counterclockwise circulation and infinity angle
$360^{\circ }/N-\unicode[STIX]{x1D719}_{\infty }$
. For example, the jets in figure 8(e,f) with
$\unicode[STIX]{x1D719}_{\infty }=48^{\circ }$
and
$\unicode[STIX]{x1D719}_{\infty }=77^{\circ }$
can be rotated
$90^{\circ }$
to yield solutions for the clockwise upper-sheet-circulation case for
$\unicode[STIX]{x1D719}_{\infty }=42^{\circ }$
and
$\unicode[STIX]{x1D719}_{\infty }=13^{\circ }$
.
3.2.6 Jetting and near-limit angles
As
$\unicode[STIX]{x1D719}_{\infty }$
approaches the limit angle
$\unicode[STIX]{x1D719}_{\infty }^{\lim }$
from below the cusp attains the Y shape in figure 8(b); above the limit the sheet and mirror image separate, forming a spiral-ended jet (figure 8
c–f).
A natural cusp–jet transition criterion is denominator
$\unicode[STIX]{x1D707}-1+e_{1}$
of the cusp exponent
$\unicode[STIX]{x1D6FC}$
in (2.45) crossing zero. But this is easily refuted: with
$N=3$
symmetry
$e_{1}=0$
so that transition would always occur at
$\unicode[STIX]{x1D707}=1$
regardless of
$\unicode[STIX]{x1D719}_{\infty }$
, clearly at odds with numerical observations (figure 10). Transition is a non-local process, influenced strongly by the ambient flow especially near the Y point rather than the cusp. Better criteria may be based on conservation of mass arguments;
$\boldsymbol{q}$
near the Y becomes very small as
$\unicode[STIX]{x1D719}_{\infty }{\nearrow}\unicode[STIX]{x1D719}_{\infty }^{\lim }$
, but we were not able to discern numerically whether it reaches zero at transition.
The value of
$\unicode[STIX]{x1D719}_{\infty }^{\lim }$
for
$\unicode[STIX]{x1D707}=1.3$
is near
$25.5^{\circ }$
, but already at
$\unicode[STIX]{x1D719}_{\infty }=48^{\circ }$
(figure 8
e), and even more so at
$\unicode[STIX]{x1D719}_{\infty }=35^{\circ }$
(figure 8
c,d) which is still almost
$10^{\circ }$
above the limit angle, the jets are very small. Distance between the spiral centres is roughly proportional to
$\unicode[STIX]{x1D719}_{\infty }-\unicode[STIX]{x1D719}_{\infty }^{\lim }$
squared. This makes computing jet solutions near the limit angle even more challenging than computing vortex cusps.
The flaws of the piecewise linear approximation occur in any flow with close and nearly parallel sheet and mirror images (even if they do not terminate in a joint cusp). In particular cusp–jet transition is predicted incorrectly, even showing coexistence in some
$\unicode[STIX]{x1D719}_{\infty }$
regions and hence spurious non-uniqueness; accurate calculations require at least quadratic degree.
3.3 Vortex cusps in Mach reflection
Henderson et al. (Reference Henderson, Vasilev, Ben-Dor and Elperin2003) and Vasilev et al. (Reference Vasilev, Ben-Dor, Elperin and Henderson2004) have begun a theoretical and numerical study on jetting in Mach reflection of shock waves for compressible flow, whose quasi-incompressible regions correspond to
$\unicode[STIX]{x1D707}=1$
in our context. It is well established (Henderson & Menikoff Reference Henderson and Menikoff1998; Elling Reference Elling2019a
) that the vortex sheets generated by single Mach reflection (figure 3
c) have counterclockwise circulation. This is very important because it allows, according to our analysis, for the presence of cusps even if the ambient flow generates only weak strain coefficients
$e_{1}$
.
Treating sheet–wall interaction as essentially incompressible allows vortex sheet numerics for better resolution. We found that when jets occur they may be small enough to be obscured by numerical or physical dissipation. The problem is compounded by the rather large cusp exponents that result from near zeros of the denominator
$\unicode[STIX]{x1D707}-1+e_{1}$
in (2.45) when
$\unicode[STIX]{x1D707}=1$
and
$e_{1}>0$
is small; such exponents cause most of the cusp to be very close to the wall, making it appear at large scale that the sheet meets the wall at a positive angle. Many flows that appear to be cusps in shock capturing or with boundary layers present are likely jets in disguise, especially if the outer part of the sheet forms a large angle
$\unicode[STIX]{x1D719}_{\infty }$
with the wall and the strain coefficient
$e_{1}$
is small and unable to compensate.
Acknowledgement
This research was partially supported by Taiwan MOST grant 105-2115-M-001-007-MY3.
