1. Introduction
In this paper we study the steady-state nonlinear flow of an ideal fluid past a submerged line vortex. As the vortices have fixed depth and horizontal displacement, they reduce to point vortices in the two-dimensional flows considered. The inviscid and incompressible fluid of infinite depth is assumed to be irrotational everywhere, except at the point vortices themselves. For a flow in the complex 
$z = x + \mathrm {i} y$ plane, with a vortex at 
$z = z_1$, the complex potential behaves as
\begin{equation} f = \phi + \mathrm{i} \psi \sim c z -\frac{\mathrm{i} \varGamma}{2 {\rm \pi}} \log(z-z_1), \end{equation}
where 
$\varGamma$ is the circulation of the vortex, and the background flow is of speed 
$c$. The non-dimensional system is then characterised by two key parameters: 
$\varGamma _{{c}}=\varGamma /(cH)$, relating vortex strength(s) to inertial effects, and the Froude number, 
$F=c/\sqrt {g H}$, relating inertial effects to gravitational effects. Here, 
$H$ is the depth of the point vortex and 
$g$ is the constant acceleration due to gravity.
The study of such vortex-driven potential flows is complicated by the following fact. The solution of two-dimensional ideal fluid-flow problems involves finding the velocity potential, 
$\phi$, and streamfunction, 
$\psi$, in terms of the coordinates 
$x$ and 
$y$, in the functional form of 
$f(z)$. However, it is often convenient to invert this dependence, instead calculating 
$z(f)$, so that the physical variables now have the forms 
$x(\phi,\psi )$ and 
$y(\phi,\psi )$. In this formulation, the free surface is a streamline along which 
$\psi$ is constant, so that the free surface is parametrised by 
$x(\phi )$ and 
$y(\phi )$. However, near the point vortex, the local behaviour (1.1) cannot be inverted analytically to give 
$z(f)$. This motivated the work of Forbes (Reference Forbes1985), who reformulated the boundary-integral formulation in terms of a free-surface arclength, 
$s$, and a more complex set of governing equations results.
The imposition of a uniform stream as 
$x \to -\infty$ results in the generation of downstream free-surface waves, as shown in figure 1(a). As hinted at in the preliminary numerical investigations of Forbes (Reference Forbes1985), the wave amplitude tends to zero as 
$F \to 0$. In this work, we confirm this behaviour and demonstrate, both numerically and analytically, that the amplitude is exponentially small in the low-Froude-number limit. For instance, the amplitude versus 
$1/F^2$ graph shown in figure 2 demonstrates the fit between our asymptotic predictions of § 3 and the numerical results of § 4. We note that this theory is nonlinear in the vortex strength, 
$\varGamma _{{c}}$, and the assumption of small 
$\varGamma _{{c}}$ need not apply.
Figure 1. The two physical regimes of underlying point vortices considered in this paper are shown. (a) A single point vortex with circulation 
$\varGamma$ is placed within the fluid. (b) Two point vortices, each with circulation 
$\varGamma$, are located at the same depth within the fluid. These solutions have been computed using the numerical scheme detailed in § 4.
Figure 2. The amplitude, 
$\bar {y}$, of the free-surface waves is shown for 
$\log (\bar {y})$ versus 
$1/F^2$ for the analytical (line) and numerical (dots) solutions of § 3 and § 4. These have a fixed value of the non-dimensional vortex strength, 
$\varGamma _{{c}}=0.25$. The graph confirms the exponential smallness of the waves. The solid line has a gradient of 
$\approx 0.7395$, computed using the exponential asymptotics theory of § 4.
The purpose of this paper is thus to characterise the formation of water waves using the framework of exponential asymptotics. We show that these exponentially small waves smoothly switch on as the fluid passes beyond the vortex, resulting in oscillations as 
$x \to \infty$ in the far field. When two submerged vortices are considered, the waves switched on due to each of the vortices may be out of phase with one another and cancel for certain values of the Froude number. This yields trapped waves between the vortices, and a free surface whose derivative decays to zero as 
$x \to \infty$. A trapped wave solution is depicted in figure 1(b). This phenomenon of trapped waves has previously been studied both for obstructions within the fluid and for flows of finite depth past lower topography. For instance, both Gazdar (Reference Gazdar1973) and Vanden-Broeck & Tuck (Reference Vanden-Broeck and Tuck1985) detected these numerically for flows over a specified lower topography. More recent works, such as those by Dias & Vanden-Broeck (Reference Dias and Vanden-Broeck2004), Hocking, Holmes & Forbes (Reference Hocking, Holmes and Forbes2013) and Holmes et al. (Reference Holmes, Hocking, Forbes and Baillard2013), have focused on detecting parameter values for which these trapped wave solutions occur in various formulations.
The work in this paper provides a first step towards extending many of the existing ideas and techniques of exponential asymptotics, previously developed for purely gravity- or capillary-driven waves (e.g. Chapman & Vanden-Broeck Reference Chapman and Vanden-Broeck2002, Reference Chapman and Vanden-Broeck2006) to wave phenomena with vortices. As noted above, because the governing equations require an alternative formulation to that originally developed by Miksis, Vanden-Broeck & Keller (Reference Miksis, Vanden-Broeck and Keller1981), the asymptotic formulation we present here can be extended to other wave–structure interactions where the more general arclength formulation of the water-wave equations is required. In addition, there has been significant recent interest in the study of water-wave phenomena with dominant vorticity effects, and we reference the recent extensive survey by Haziot et al. (Reference Haziot, Hur, Strauss, Toland, Wahlén, Walsh and Wheeler2022) and references therein. The exponential asymptotics techniques developed in this work can also be extended to situations where capillary ripples are forced on the surface of steep vortex-driven waves. The leading-order solution in these asymptotic regimes would then be known analytically from the works of Crowdy & Nelson (Reference Crowdy and Nelson2010), Crowdy & Roenby (Reference Crowdy and Roenby2014) and Crowdy (Reference Crowdy2023), for example. We shall discuss these and other exciting future directions in § 6.
2. Mathematical formulation and outline
We consider the typical configurations shown in figure 1. Following Forbes (Reference Forbes1985), in non-dimensional form, the system is formulated in terms of the arclength, 
$s$, along the free surface, with unknown velocity potential, 
$\phi = \phi (s)$, and free-surface positions, 
$(x(s), y(s))$. Then, the governing equations are given by Bernoulli's equation, an arclength relation between 
$x$ and 
$y$, and a boundary-integral equation.
For a single submerged point vortex at 
$(x,y)=(0,-1)$, the three equations are
\begin{gather} \frac{F^2}{2}[\phi^{\prime}(s)]^2+y(s)=\frac{F^2}{2}, \end{gather}
\begin{gather}{}[x^{\prime}(s)]^2+[y^{\prime}(s)]^2=1, \end{gather}
\begin{gather}\phi^{\prime}(s) x^{\prime}(s) - 1 = \frac{\varGamma_{{c}}}{\rm \pi} \frac{y(s)+1}{[x(s)]^2 + [y(s)+1]^2} + \mathcal{I}[x,y,\phi^{\prime}]. \end{gather}
In the above, two non-dimensional parameters appear: the Froude number, 
$F$, and the vortex strength, 
$\varGamma _{{c}}$, defined by
\begin{equation} F= \frac{c}{\sqrt{g H}} \quad \text{and} \quad \varGamma_{{c}}=\frac{\varGamma}{cH}. \end{equation}
Here, 
$c$ is the speed of the fluid, 
$H$ is the depth of the submerged point vortex, 
$g$ is the constant acceleration due to gravity and 
$\varGamma$ is the circulation of the point vortex. Furthermore, we have also introduced 
$\mathcal {I}$ as the nonlinear principal-value integral defined by
\begin{equation} \mathcal{I}[x,y,\phi^{\prime}]=\frac{1}{\rm \pi}\,{\int\hskip -1,05em -\,}_{-\infty}^{\infty} \frac{[\phi^{\prime}(t)-x^{\prime}(t)][y(t)-y(s)]+y^{\prime}(t)[x(t)-x(s)]} {[x(t)-x(s)]^2+[y(t)-y(s)]^2} \, \mathrm{d} t. \end{equation}When the configuration with two point vortices is considered in § 3.5, the boundary-integral equation (2.1c) will need to be modified to (3.23).
2.1. Analytic continuation
In the exponential asymptotics procedure of § 3, we study the exponentially small terms that display the Stokes phenomenon across Stokes lines of the problem. These Stokes lines originate from singularities of the leading-order asymptotic solution, which are located in the analytic continuation of the domain, the arclength 
$s$. The analytic continuation of the governing equations (2.1a)–(2.1c) is studied in this section.
We now analytically continue the domain 
$s \mapsto \sigma$, where 
$\sigma \in \mathbb {C}$. Bernoulli's equation (2.1a) and the arclength relation (2.1b) may be analytically continued in a straightforward manner, with all dependence on 
$s$ replaced by the complex-valued variable 
$\sigma$. The analytic continuation of the boundary-integral equation (2.1c) is more complicated, due to the principal-value integral 
$\mathcal {I}$ defined in (2.3). The analytic continuation of this integral is given by
\begin{equation} \mathcal{I}[x,y,\phi^{\prime}]= \hat{\mathcal{I}}[x,y,\phi^{\prime}] - a \mathrm{i} \phi^{\prime}(\sigma)y^{\prime}(\sigma), \end{equation}
where 
$a=\pm 1$ denotes the direction of analytic continuation into 
$\text {Im}[\sigma ]>0$ or 
$\text {Im}[\sigma ]<0$, respectively, and 
$\hat {\mathcal {I}}$ is the complex-valued integral. Equation (2.4) may be verified by taking the limit of either 
$\text {Im}[\sigma ] \to 0^+$ or 
$\text {Im}[\sigma ] \to 0^-$, which yields half of a residue contribution associated with the singular point at 
$t=s$ of the integrand.
Substitution of (2.4) into (2.1c) then yields the analytically continued equations, given by
\begin{gather} \frac{F^2}{2}[\phi^{\prime}(\sigma)]^2+y(\sigma)=\frac{F^2}{2}, \end{gather}
\begin{gather}{}[x^{\prime}(\sigma)]^2+[y^{\prime}(\sigma)]^2=1, \end{gather}
\begin{gather}\phi^{\prime}(\sigma) x^{\prime}(\sigma) -1 + a \mathrm{i} \phi^{\prime}(\sigma)y^{\prime}(\sigma) = \frac{\varGamma_{{c}}}{\rm \pi} \frac{y(\sigma)+1}{[x(\sigma)]^2 + [y(\sigma)+1]^2}+ \hat{\mathcal{I}}[x,y,\phi^{\prime}]. \end{gather}The analytic continuation for situations with multiple point vortices is done similarly, with the only difference being the inclusion of additional point vortices in (2.5c).
2.2. Outline of the paper
In this work, we will consider the following two regimes depicted in figure 1:
(a) A single submerged point vortex, which is the formulation originally considered by Forbes (Reference Forbes1985). Imposing free-stream conditions as
$x \to -\infty$ results in surface waves generated by the vortex. Their amplitude is exponentially small as $F \to 0$. This is the limit considered by Chapman & Vanden-Broeck (Reference Chapman and Vanden-Broeck2006) in the absence of vortical effects.(b) Two submerged point vortices of the same circulation. For certain critical values of the Froude number,
$F$, the resultant waves are confined to lie between the two vortices. The amplitude of these is also exponentially small as $F \to 0$.
We begin in § 3 by determining these exponentially small waves using the techniques of exponential asymptotics. This relies on the optimal truncation of an algebraic asymptotic series for small Froude number, 
$F$, and deriving the connection of this to the Stokes phenomenon that acts on the exponentially small waves. The case for two submerged point vortices is then studied in § 3.5, where we derive the critical values of the Froude number for which the waves are trapped. Numerical solutions are computed in § 4, where comparison occurs with the exponential asymptotics predictions for the single-vortex and double-vortex cases.
3. Exponential asymptotics
3.1. Early orders of the solution
We begin by considering the following asymptotic expansions, in powers of 
$F^2$, for the solutions, which are given by
\begin{equation} x(\sigma)= \sum_{n=0}^{\infty} F^{2n}x_n(\sigma), \quad y(\sigma)= \sum_{n=0}^{\infty} F^{2n}y_n(\sigma), \quad \phi^{\prime}(\sigma)= \sum_{n=0}^{\infty}F^{2n}\phi_n^{\prime}(\sigma). \end{equation}
Substitution of expansions (3.1a–c) into (2.5a)–(2.5c) yields at leading order three equations for the unknowns 
$x_0$, 
$y_0$ and 
$\phi _0^{\prime }$. The first of these, Bernoulli's equation (2.5a), yields 
$y_0(\sigma )=0$. This may be substituted into the second equation, (2.5b), to find 
$(x_0^{\prime })^2=1$, for which we consider 
$x_0^{\prime }=1$ without any loss of generality. This may be integrated to find 
$x_0=\sigma$, where the constant of integration has been chosen to set the origin at 
${x_0(0)=0}$. Next, 
$\phi _0^{\prime }$ is determined from (2.5c). Since 
$y_0=0$, the integral 
$\hat {\mathcal {I}}$ does not enter the leading-order equation. This yields the leading-order solutions as
\begin{equation} y_0(\sigma)=0, \quad x_0(\sigma)=\sigma, \quad \phi^{\prime}_0(\sigma)=1+\frac{\varGamma_{{c}}}{\rm \pi}\frac{1}{(1+\sigma^2)}. \end{equation} Note that there is a singularity in 
$\phi _0^{\prime }$ above whenever 
$\sigma ^2=-1$. This corresponds to the point vortex within the fluid at 
$\sigma =-\mathrm {i}$, as well as another singularity at 
$\sigma = \mathrm {i}$, which will produce a complex-conjugate contribution to the exponentially small solution along the free surface.
Next, at order 
$O(F^2)$, 
$y_1$ is found explicitly from (2.5a). We then find the equation 
$x_1^{\prime }=0$ from (2.5b), and 
$\phi _1$ is determined explicitly from (2.5c). This yields
\begin{equation} \left. \begin{gathered} y_1(\sigma)=\tfrac{1}{2}[1-(\phi_0^{\prime})^2], \quad x_1(\sigma)=0,\\ \phi^{\prime}_1(\sigma) ={-}a \mathrm{i} \phi_0^{\prime}y_1^{\prime} +\frac{\varGamma_{{c}}(\sigma^2-1)}{{\rm \pi}(1+\sigma^2)^2}y_1 +\hat{\mathcal{I}}_1(\sigma), \end{gathered} \right\} \end{equation}
where 
$\hat {\mathcal {I}}_1$ is the 
$O(F^2)$ component of the complex-valued integral 
$\hat {\mathcal {I}}$, originally defined along the real axis in (2.3).
3.2. Late-term divergence
Our derivation of the exponentially small terms and associated Stokes phenomenon of § 3.4 requires knowledge of the late terms of the solution expansion (3.1a–c), i.e. 
$x_n$, 
$y_n$ and 
$\phi _n^{\prime }$, as 
$n \to \infty$. We begin by determining the 
$O(F^{2n})$ components of (2.5a)–(2.5c). The late terms of Bernoulli's equation are given by
\begin{equation} y_n + \phi_0^{\prime} \phi_{n-1}^{\prime} + \phi_1^{\prime}\phi_{n-2}^{\prime} +\cdots =0, \end{equation}for the arclength relation we have
\begin{equation} x_0^{\prime}x_n^{\prime} + x_1^{\prime}x_{n-1}^{\prime}+\cdots + y_1^{\prime}y_{n-1}^{\prime}+ y_2^{\prime}y_{n-2}^{\prime}+ \cdots=0, \end{equation}and finally the boundary-integral equation yields
\begin{align} & x_0^{\prime}\phi_n^{\prime} + x_1^{\prime}\phi_{n-1}^{\prime} + \phi_0^{\prime}x_n^{\prime} + \cdots + a \mathrm{i} [\phi_0^{\prime}y_n^{\prime} + \phi_1^{\prime}y_{n-1}^{\prime} + y_1^{\prime}\phi_{n-1}^{\prime} + \cdots]\nonumber\\ &\quad + \frac{\varGamma_{{c}}}{\rm \pi} \left[ \frac{y_n}{1+x_0^2} - \frac{2y_n}{(1+x_0^2)^2}+\cdots \right] -\hat{\mathcal{I}}_n(\sigma)=0. \end{align}
In (3.4a)–(3.4c) above, only the terms that will appear at the first two orders of 
$n$ as 
$n \to \infty$ have been included.
In (3.4c), the 
$O(F^{2n})$ component of the complex-valued integral, 
$\hat {\mathcal {I}}$, has been denoted by 
$\hat {\mathcal {I}}_n$. The dominant components of this integral, as 
$n \to \infty$, require the integration of late-term asymptotic solutions that are either a function of the real-valued integration domain, such as 
$y_n(t)$, or a function of the complex domain, such as 
$y_n(\sigma )$. The first of these, 
$y_n(t)$, is integrated along the real-valued free surface, away from any singular behaviour. It is thus subdominant to the other terms appearing in (3.4c). This is analogous to neglecting the late terms of the complex-valued Hilbert transform in similar free-surface problems in exponential asymptotics (cf. Xie & Tanveer Reference Xie and Tanveer2002; Chapman & Vanden-Broeck Reference Chapman and Vanden-Broeck2002, Reference Chapman and Vanden-Broeck2006). All that remains is to integrate the components of 
$\hat {\mathcal {I}}_n$ that involve late-term solutions evaluated in the complex-valued domain. Of these, only that involving 
$y_n(\sigma )$ appears in the two leading orders, as 
$n \to \infty$, of (3.4c). This component is given by
\begin{equation} \hat{\mathcal{I}}_n \sim{-}\frac{y_n(\sigma)}{\rm \pi}\int_{-\infty}^{\infty} \frac{\phi_0^{\prime}(t)-1}{(t-\sigma)^2}\,\mathrm{d}t ={-}\frac{\varGamma_{{c}}}{\rm \pi}\frac{y_n(\sigma)}{(\sigma+a\mathrm{i})^2}, \end{equation}
for which the integral was evaluated by substituting for 
$\phi _0^{\prime }$ from (3.2a–c). Note that integration of 
$y_n(\sigma )$ was not required due to the lack of any dependence on the domain of integration, 
$t$.
Recall that the leading-order solutions were singular at 
$\sigma = \pm \mathrm {i}$. For each of the three solution expansions, this singularity first appeared in 
$\phi ^{\prime }_0$, 
$y_1$ and 
$x_2$. Since successive terms in the asymptotic expansion involve differentiation of previous terms (for instance, (3.4a) for 
$y_n$ involves 
$\phi ^{\prime }_{n-1}$, whose determination in (3.4c) requires knowledge of 
$y_{n-1}^{\prime }$), the strength of this singularity will grow as we proceed into the asymptotic series. Furthermore, this growing singular behaviour will also lead to the divergence of the late-term solutions as 
$n \to \infty$, which we capture analytically with the factorial-over-power ansatzes of
\begin{equation} x_n \sim X(\sigma)\frac{\varGamma(n+\alpha-1)}{[\chi(\sigma)]^{n+\alpha-1}}, \quad y_n \sim Y(\sigma)\frac{\varGamma(n+\alpha)}{[\chi(\sigma)]^{n+\alpha}}, \quad \phi_n \sim \varPhi(\sigma)\frac{\varGamma(n+\alpha)} {[\chi(\sigma)]^{n+\alpha}}. \end{equation}
Here, 
$\alpha$ is a constant, 
$\chi$ is the singulant function that will capture the singular behaviour of the solution at 
$\sigma =\pm \mathrm {i}$, and 
$X$, 
$Y$ and 
$\varPhi$ are functional prefactors of the divergent solutions. It can be seen from the dominant balance as 
$n \to \infty$ of (3.4a) and (3.4b) that 
$x_{n+1} =O( y_{n})$ and 
$y_n=O( \phi _{n})$, which has motivated our precise ordering in 
$n$ in the ansatzes (3.6a–c).
Substitution of ansatzes (3.6a–c) into the 
$O(F^{2n})$ equations (3.4a)–(3.4c) yields, at leading order in 
$n$, the three equations
\begin{equation} Y- \phi_0^{\prime}\chi^{\prime}\varPhi=0, \quad \chi^{\prime}(X+ y_1^{\prime} Y)=0, \quad \chi^{\prime}(\varPhi+a \mathrm{i} \phi_0^{\prime}Y)=0. \end{equation}
While the last two of these equations permit the solution 
$\chi ^{\prime }=0$, this is unable to satisfy the first equation in (3.7a–c). The remaining solutions can be solved to give 
$\chi ^{\prime }= a \mathrm {i} (\phi _0^{\prime })^{-2}$, which we integrate to find
\begin{equation} \chi_a(\sigma)=a \mathrm{i} \int_{a \mathrm{i}}^{\sigma} \left[1+\frac{\varGamma_{{c}}}{\rm \pi}\frac{1}{(1+t^2)}\right]^{{-}2}\mathrm{d}t. \end{equation}
Here, we have introduced the notation 
$\chi _a=\chi$, where 
$a=\pm 1$, to discern between each singulant generated by the two singular points of 
$\phi _0^{\prime }$, which are given by 
$\sigma =\mathrm {i}$ and 
$\sigma = -\mathrm {i}$. The starting point of integration in (3.8) is 
$\sigma = \pm \mathrm {i}$ to ensure that 
$\chi _a(a \mathrm {i})=0$. This condition is required in order to match with an inner solution near this singular point. Integration of (3.8) yields
\begin{align} \chi_a(\sigma) & = a \mathrm{i} \left[\sigma + \frac{\varGamma_{{c}}^2\sigma}{2(\varGamma_{{c}}+{\rm \pi})({\rm \pi} \sigma^2 + \varGamma_{{c}} + {\rm \pi})}-\frac{\varGamma_{{c}} (3 \varGamma_{{c}}+4 {\rm \pi})} {2\sqrt{\rm \pi}(\varGamma_{{c}}+{\rm \pi})^{3/2}} \tan^{{-}1}{\left(\frac{ \sqrt{\rm \pi} \sigma}{\sqrt{\varGamma_{{c}}+{\rm \pi}}}\right)} \right]\nonumber\\ & \quad +1+\frac{\varGamma_{{c}}}{2(\varGamma_{{c}}+{\rm \pi})} -\frac{\varGamma_{{c}}(3 \varGamma_{{c}}+4{\rm \pi})} {2 \sqrt{\rm \pi}(\varGamma_{{c}}+{\rm \pi})^{3/2}} \tanh^{{-}1}{\left( \frac{\sqrt{\rm \pi}}{\sqrt{\varGamma_{{c}}+{\rm \pi}}}\right)}. \end{align}3.3. Solution of the late-term amplitude equations
We now determine the amplitude functions, 
$\varPhi$, 
$X$ and 
$Y$, of the late-term solutions. Note that, if one of these amplitude functions is known, then the other two may be determined by the last two equations in (3.7a–c). Thus, only one equation is required for the amplitude functions, which we find at the next order of 
$n$ in the late-term equation (3.4a). This equation is given by
\begin{equation} \phi_0^{\prime} \varPhi^{\prime} = \phi_1^{\prime} \chi^{\prime} \varPhi, \end{equation}which may be integrated to find the solution
\begin{equation} \varPhi(\sigma)=\varLambda \exp\left( a \mathrm{i} \int_{0}^{\sigma} \frac{\phi_1^{\prime}(t)}{[\phi_0^{\prime}(t)]^3}\,\mathrm{d}t\right). \end{equation}
In the above, 
$\varLambda$ is a constant of integration, which is determined by matching with an inner solution near the singular points 
$\sigma =a\mathrm {i}$. Once 
$\varPhi$ is known, the remaining amplitude functions are determined by the equations 
$Y= a \mathrm {i} (\phi _0^{\prime })^{-1} \varPhi$ and 
$X= a \mathrm {i} \phi _0^{\prime \prime } \varPhi$.
We now calculate the constant, 
$\alpha$, that appears in the factorial-over-power ansatzes (3.6a–c). This is determined by ensuring that the singular behaviour, as 
$\sigma \to a \mathrm {i}$, of each ansatz is consistent with the anticipated singular behaviours of
\begin{equation} x_n = O( (\sigma - a\mathrm{i})^{1-3n} ), \quad y_n = O( (\sigma - a\mathrm{i})^{1-3n} ), \quad \phi_n = O( (\sigma - a\mathrm{i})^{{-}3n} ). \end{equation}
In taking the inner limit of 
$\varPhi$ from (3.11), we have 
$\varPhi = O(\sigma -a \mathrm {i})^{3/2}$. Furthermore, since 
$\chi = O((\sigma - a \mathrm {i})^{3})$, derived later in (A5), equating the power of the singularities for 
$\phi _n$ between the ansatzes (3.6a–c) and the anticipated singular behaviour above in (3.12a–c) yields the value of 
$\alpha =1/2$. The constant of integration, 
$\varLambda$, that appears in solution (3.11) for the amplitude function, 
$\varPhi$, is derived in Appendix A by matching the inner limit of the divergent solution, 
$\phi _n$, with an inner solution at 
$\sigma =a \mathrm {i}$. This yields
\begin{equation} \alpha = \frac{1}{2} \quad \text{and} \quad \varLambda ={-}\frac{a \mathrm{i} \mathrm{e}^{-\mathcal{P}(a \mathrm{i})}}{3 \sqrt{3}} \lim_{n \to \infty} \left(\frac{\hat{\phi}_n}{\varGamma(n+\alpha+1)}\right), \end{equation}
where 
$\hat {\phi }_n$, determined via recurrence relation (A10), is a constant appearing in the series expansion for the outer limit of the inner solution for 
$\phi$, and 
$\mathcal {P}(\sigma )$ is defined in (A13).
To conclude, the late terms of the asymptotic expansions (3.1a–c) diverge in a factorial-over-power manner specified by the ansatzes (3.6a–c). Evaluation of this divergence requires the constants 
$\alpha$ and 
$\varLambda$ from (3.13a,b), as well as the singulant function 
$\chi (\sigma )$ from (3.9) and the amplitude function 
$\varPhi (\sigma )$ from (3.11). These will be required in the derivation of the exponentially small terms considered in the next section.
3.4. Stokes smoothing and Stokes lines
The exponentially small components of the solutions are now determined. We truncate the asymptotic expansions (3.1a–c) at 
$n=N-1$ and consider the remainders, yielding
\begin{equation} x= \underbrace{\sum_{n=0}^{N-1} F^{2n} x_n}_{x_r} + \bar{x}, \quad y= \underbrace{\sum_{n=0}^{N-1} F^{2n} y_n}_{y_r} + \bar{y}, \quad \phi^{\prime}=\underbrace{ \sum_{n=0}^{N-1}F^{2n}\phi_n^{\prime}}_{\phi^{\prime}_r}+ \bar{\phi}, \end{equation}
where the truncated asymptotic expansions have been denoted by 
$x_r$, 
$y_r$ and 
$\phi _r^{\prime }$. When 
$N$ is chosen optimally at the point at which the divergent expansions reorder as 
$n \to \infty$, given by
\begin{equation} N \sim \frac{\lvert \chi \rvert}{F^2} + \rho, \end{equation}
where 
$0\leq \rho <1$ to ensure that 
$N$ is an integer, the remainders to the asymptotic expansions (3.14a–c) will be exponentially small.
Equations for these remainders are found by substituting the truncated expansions (3.14a–c) into the analytically continued equations (2.5a)–(2.5c). These are given by
\begin{gather} (F^2\phi_0^{\prime}+F^4\phi_1^{\prime})\bar{\phi}^{\prime}+\bar{y}={-}\xi_{{a}}, \end{gather}
\begin{gather}2 \bar{x}^{\prime} +2F^2 y_1^{\prime}\bar{y}^{\prime} ={-}\xi_{{b}}, \end{gather}
\begin{gather}\bar{\phi}^{\prime}+a \mathrm{i} \phi_0^{\prime} \bar{y}^{\prime}={-}\xi_{{c}}. \end{gather}
In equations (3.16) above, nonlinear terms such as 
$\bar {x}^2$ were neglected, as they will be exponentially subdominant. In anticipating that 
$\bar {x}=O(F^2\bar {y})=O(F^2\bar {\phi })$, terms of the first two orders of 
$F^2$ have been retained on the left-hand side of (3.16a). Motivated by the late-term analysis, in which equations for the amplitude functions were obtained at leading order for the last two governing equations, we have only retained the leading-order terms in (3.16b) and (3.16c). Furthermore, the forcing terms introduced in (3.16) are defined by
\begin{equation} \left. \begin{gathered} \xi_{{a}}=\frac{F^2}{2}(\phi_r^{\prime})^2+y_r-\frac{F^2}{2}, \quad \xi_{{b}}=(x_r^{\prime})^2+(y_r^{\prime})^2-1,\\ \xi_{{c}}=\phi_r^{\prime}x_r^{\prime}-1+a \mathrm{i}\phi_r^{\prime}y_r^{\prime} -\frac{\varGamma_{{c}}}{\rm \pi}\frac{y_r+1}{(x_r)^2+(y_r+1)^2} -\hat{\mathcal{I}}[x_r,y_r,\phi_r^{\prime}]. \end{gathered} \right\} \end{equation}
Since each order of these forcing terms will be identically zero up to and including 
$O(F^{2(N-1)})$, each of equations (3.17) will be of 
$O(F^{2N})$. Only knowledge of 
$\xi _{{a}}$ will be required in the Stokes smoothing procedure of this section, and the leading-order component is given by
\begin{equation} \xi_{{a}}\sim \phi_0^{\prime}\phi_{N-1}^{\prime}F^{2N}. \end{equation} Homogeneous solutions to (3.16), for which the forcing terms on the right-hand sides are omitted, are 
$\bar {x} \sim F^2 X\mathrm {e}^{-\chi /F^2}$, 
$\bar {y} \sim Y\mathrm {e}^{-\chi /F^2}$ and 
$\bar {\phi } \sim \varPhi \mathrm {e}^{-\chi /F^2}$, where the singulant 
$\chi$ and amplitude functions 
$X$, 
$Y$ and 
$\varPhi$ satisfy the same equations as those found for the late-term solutions in § 3.2. Next, we solve for the particular solutions of (3.16) through variation of parameters by multiplying the homogeneous solutions by an unknown function, 
$\mathcal {S}(\sigma )$, giving
\begin{equation} \left. \begin{gathered} \bar{x} \sim \mathcal{S}(\sigma)F^2X(\sigma)\mathrm{e}^{-\chi(\sigma)/F^2}, \\ \bar{y} \sim \mathcal{S}(\sigma) Y(\sigma)\mathrm{e}^{-\chi(\sigma)/F^2}, \\ \bar{\phi} \sim \mathcal{S}(\sigma)\varPhi(\sigma)\mathrm{e}^{-\chi(\sigma)/F^2}, \end{gathered} \right\} \end{equation}
where 
$Y=a \mathrm {i} \varPhi / \phi _0^{\prime }$ and 
$X=-y_1^{\prime }Y$. The function 
$\mathcal {S}$ is called the Stokes multiplier, as it will display the Stokes phenomenon across Stokes lines of the problem, which is demonstrated next.
An equation for 
$\mathcal {S}$ is obtained by substituting (3.19) into (3.16a), yielding 
$F^2 \phi _0^{\prime }\varPhi \mathrm {e}^{-\chi /F^2}\mathcal {S}^{\prime }(\sigma ) \sim - \xi _{{a}}$. In substituting for the dominant behaviour of 
$\xi _{{a}}$ from (3.18) and the factorial-over-power divergence of 
$\phi ^{\prime }_{N-1}$ from (3.6a–c), we change derivatives of 
$\mathcal {S}$ from 
$\sigma$ to 
$\chi$ to find
\begin{equation} \frac{\mathrm{d}\mathcal{S}}{\mathrm{d} \chi} \sim \frac{\varGamma(N+\alpha)}{\chi^{N+\alpha}}F^{2(N-1)}\mathrm{e}^{\chi/F^2}. \end{equation}
In expanding as 
$N \to \infty$, and substituting for 
$N \sim \lvert \chi \rvert /F^2 + \rho$ from (3.15), the right-hand side of (3.19) is seen to be exponentially small, except for in a boundary layer close to contours satisfying
\begin{equation} \text{Im}[\chi]=0 \quad \text{and} \quad \text{Re}[\chi]>0. \end{equation}
These are the Stokes line conditions originally derived by Dingle (Reference Dingle1973). Across the Stokes lines, the solution for the Stokes multiplier 
$\mathcal {S}$,
\begin{equation} \mathcal{S}(\sigma)=S_a+\frac{\sqrt{2 {\rm \pi}} \mathrm{i}}{F^{2 \alpha}} \int_{-\infty}^{\sqrt{\lvert\chi\rvert}({\arg(\chi)}/{F})} \exp({-}t^2/2) \,\mathrm{d}t, \end{equation}
rapidly varies from the constant 
$S_a$ to 
$S_a + 2 {\rm \pi}\mathrm {i} / F^{2\alpha }$. This is the Stokes phenomenon, and the contours satisfying the Dingle conditions (3.21a,b) are shown in figure 3 to lie along the imaginary axis. For the one-vortex case studied in this section, the upstream condition as 
$\textrm {Re}[\sigma ] \to -\infty$ requires that 
$S_1=0$ and 
$S_{-1}=-2 {\rm \pi}\mathrm {i} /F^{2 \alpha }$.
Figure 3. The Stokes lines (bold) lie along the imaginary axis between the two singular points of 
$\sigma =-\mathrm {i}$ and 
$\sigma =\mathrm {i}$. Branch cuts are shown with a wavy line.
3.5. Trapped waves generated by two submerged vortices
We have so far studied the case of a single submerged point vortex. When multiple point vortices are placed within the fluid, the only change is to the boundary-integral equation, previously specified in (2.5c) for a single vortex. In this section we study the formulation of two submerged point vortices of the same non-dimensional strength, 
$\varGamma _{{c}}$, located at 
$z=x+\mathrm {i} y = \pm \lambda - \mathrm {i}$, for which the analytically continued boundary-integral equation is given by
\begin{align} & \phi^{\prime}(\sigma)x^{\prime}(\sigma) - 1 + a \mathrm{i} \phi^{\prime}(\sigma)y^{\prime}(\sigma) \nonumber\\ & \quad = \frac{\varGamma_{{c}}}{\rm \pi} \left[\frac{ y(\sigma)+1}{[x(\sigma)-\lambda]^2 + [y(\sigma)+1]^2} + \frac{ y(\sigma)+1}{[x(\sigma)+\lambda]^2 + [y(\sigma)+1]^2}\right] + \hat{\mathcal{I}}[x,y,\phi]. \end{align} Unlike the case for a single submerged point vortex that produces waves in the far field for 
$x \to \infty$, two identical point vortices can produce solutions for which the waves are confined to lie between the vortices, 
$-\lambda <\textrm {Re}[\sigma ]<\lambda$. This occurs for critical values of the Froude number, which we now predict using the techniques of exponential asymptotics developed in the previous sections.
The first two orders of the asymptotic solution for 
$\phi$ are now given by
\begin{gather} \phi_0^{\prime}(\sigma)=1+\frac{\varGamma_{{c}}}{\rm \pi} \left[\frac{1}{1+(\sigma+\lambda)^2}+\frac{1}{1+(\sigma-\lambda)^2}\right], \end{gather}
\begin{gather}\phi_1^{\prime}(\sigma)={-} a \mathrm{i} \phi_0^{\prime}y_1^{\prime}+\frac{\varGamma_{{c}}y_1}{\rm \pi} \left[\frac{(\sigma+\lambda)^2-1}{[1+(\sigma+\lambda)^2]^2} +\frac{(\sigma-\lambda)^2-1}{[1+(\sigma-\lambda)^2]^2}\right] +\hat{\mathcal{I}}_n(\sigma), \end{gather}
which are singular at the four locations 
$\sigma =-\lambda +a \mathrm {i}$ (from the vortex at 
$z=-\lambda -\mathrm {i}$) and 
$\sigma =\lambda + a \mathrm {i}$ (from the vortex at 
$z=\lambda -\mathrm {i}$). Note that we have again defined 
$a=\pm 1$ to indicate whether 
$\textrm {Im}[\sigma ]>0$ or 
$\textrm {Im}[\sigma ]<0$. These four singular points each have associated Stokes lines, shown in figure 4. In general, the waves switched on across the first Stokes lines, emanating from the points 
$\sigma =-\lambda + a \mathrm {i}$, will be out of phase with the waves switched on across the second Stokes lines, from 
$\sigma =\lambda +a \mathrm {i}$. However, for certain values of 
$F$, the wave switched on across the first Stokes line is then switched off by the second Stokes line, yielding solutions with no waves for 
$\textrm {Re}[\sigma ]> \lambda$. An example of this trapped solution was shown earlier in figure 1(b).
Figure 4. The Stokes lines (bold) generated by the four singular points are shown.
Thus, in using the Stokes switching prediction for 
$\bar {\phi }$ shown in figure 4 and writing 
$\bar {y}= a \mathrm {i} \bar {\phi } / \phi _0^{\prime }$, we require the two contributions of
\begin{equation} \left. \begin{gathered} \bar{y}_1 \sim{-}\frac{2 {\rm \pi}}{F^{2 \alpha}\phi^{\prime}_0 }\varPhi_1(\sigma) \exp\left(-\frac{\chi_1(\sigma)}{F^2}\right) + {\rm c.c.},\\ \bar{y}_2 \sim{-}\frac{2 {\rm \pi}}{F^{2 \alpha}\phi^{\prime}_0 }\varPhi_2(\sigma) \exp\left(-\frac{\chi_2(\sigma)}{F^2}\right) + {\rm c.c.}, \end{gathered} \right\} \end{equation}
to cancel with one another for 
$\textrm {Re}[\sigma ]>\lambda$. Here, we have denoted 
$\chi _1$ and 
$\varPhi _1$ as the singulant and amplitude function arising from the 
$\sigma =-\lambda +a \mathrm {i}$ singularities, and 
$\chi _2$ and 
$\varPhi _2$ as those arising from the 
$\sigma =\lambda +a \mathrm {i}$ singularities. The first of equations (3.25), 
$\bar {y}_1$, is the contribution switched on as we pass from left to right across the Stokes lines associated with the singular points 
$\sigma =-\lambda +a\mathrm {i}$. The second, 
$\bar {y}_2$, is the contribution switched on from left to right by the Stokes lines associated with the 
$\sigma =\lambda +a\mathrm {i}$ singular point. Note that the specified contributions in (3.25) are from the 
$a=1$ contribution, and the unspecified complex-conjugate components are from that with 
$a=-1$.
We now simplify each of the expressions given in (3.25) by substituting for the amplitude functions 
$\varPhi _1$ and 
$\varPhi _2$, which satisfy the same equation as that found previously in (3.10). The only difference will be the constants of integration, which we denote by 
$\varLambda _1$ and 
$\varLambda _2$. This yields
\begin{equation} \left. \begin{gathered} \bar{\phi}_1 \sim{-}\frac{4 {\rm \pi}\lvert \varLambda_1 \rvert }{F^{2 \alpha}\phi_0^{\prime}} \exp\left(-\frac{\text{Re}[\chi_1]}{F^2}\right) \cos\left(\int_{0}^{\sigma} \frac{\phi_1^{\prime}(t)}{[\phi_0^{\prime}(t)]^3}\,\mathrm{d}t +\arg[\varLambda_1]-\frac{\text{Im}[\chi_1]}{F^2}\right),\\ \bar{\phi}_2 \sim{-}\frac{4 {\rm \pi}\lvert \varLambda_2 \rvert }{F^{2 \alpha}\phi_0^{\prime}} \exp\left(-\frac{\text{Re}[\chi_2]}{F^2}\right) \cos\left(\int_{0}^{\sigma} \frac{\phi_1^{\prime}(t)}{[\phi_0^{\prime}(t)]^3}\,\mathrm{d}t +\arg[\varLambda_2]-\frac{\text{Im}[\chi_2]}{F^2}\right). \end{gathered} \right\} \end{equation}
Through integration of 
$\chi ^{\prime }=a \mathrm {i} (\phi _0^{\prime })^{-2}$ and imposing the boundary conditions 
$\chi _1(a \mathrm {i} - \lambda )=0$ and 
$\chi _2(a \mathrm {i} + \lambda )=0$, it may be verified that, along the free surface, 
$\textrm {Im}[\sigma ]=0$, we have 
$\textrm {Re}[\chi _1]=\textrm {Re}[\chi _2]$. Furthermore, we also have 
$\lvert \varLambda _1 \rvert =\lvert \varLambda _2 \rvert$. This may be verified by matching with an inner solution, much like the procedure considered in Appendix A for the case of a single point vortex. The same leading-order inner equation emerges regardless of the number of vortices considered, and so the outer limit of the inner solution is the same as in (A11) but with functional dependence on either 
$\sigma - \lambda -a \mathrm {i}$ or 
$\sigma + \lambda - a\mathrm {i}$. Thus, the only difference encountered in the matching procedure is in the inner limit of the outer divergent solution.
The prefactors multiplying each of the cosine functions in (3.26) are identical, and the condition for them to cancel, 
$\bar {y}_1+\bar {y}_2=0$, yields
\begin{align} &\cos\left( \int_0^{\sigma} \frac{\phi_1^{\prime}(t)}{[\phi_0^{\prime}(t)]^3}\,\mathrm{d}t + \frac{\arg[\varLambda_1] + \arg[\varLambda_2]}{2} -\frac{\text{Im}[\chi_1 + \chi_2]}{2F^2} \right)\nonumber\\ &\quad\times \cos\left(\frac{\arg[\varLambda_1]-\arg[\varLambda_2]}{2} -\frac{\text{Im}[\chi_1 - \chi_2]}{2F^2} \right) =0. \end{align}
Note that, since 
$\chi _1$ and 
$\chi _2$ satisfy the same differential equation, 
$\chi ^{\prime }=a \mathrm {i} (\phi _0^{\prime })^{-2}$, originally derived in § 3.2, the only difference between them are their constants of integration. Therefore, 
$\textrm {Im}[\chi _1+\chi _2]$ will be a function of 
$\sigma$, and 
$\textrm {Im}[\chi _1-\chi _2]$ will be constant. Thus, only the second cosine component of (3.27) is capable of satisfying the identity for 
$\textrm {Re}[\sigma ]>\lambda$. Since this cosine function is zero when the argument equals 
$\pm {\rm \pi}/2$, 
$\pm 3 {\rm \pi}/2$, and so forth, we find
\begin{equation} F_k = \sqrt{\frac{\text{Im}[\chi_1 - \chi_2]}{ \arg[\varLambda_1]-\arg[\varLambda_2]+{\rm \pi}(2k+1)}}, \end{equation}
for 
$k =0,1,2,\ldots\,$, and so forth. Equation (3.28) yields the discrete values of the Froude number, 
$F_k$, for which the waves are confined to lie between the two submerged vortices.
All that remains is to evaluate 
$\textrm {Im}[\chi _1-\chi _2]$, 
$\arg [\varLambda _1]$ and 
$\arg [\varLambda _2]$. Each of these singulants is found by integrating 
$\chi ^{\prime }=a\mathrm {i} (\phi _0^{\prime })^{-2}$, where 
$\phi _0^{\prime }$ is specified in (3.24), from the corresponding singular point. We may decompose each singulant into a real-valued integral along the Stokes line, and an imaginary-valued integral along the free surface. Thus, 
$\textrm {Im}[\chi ]$ is an integral along the free surface, 
$\textrm {Im}[\sigma ]=0$, from the intersection of the Stokes line to 
$\sigma$. This yields
\begin{equation} \text{Im}[\chi_1(\sigma) - \chi_2(\sigma)] = \int_{-\lambda}^{\lambda} \left[1+\frac{\varGamma_{{c}}}{\rm \pi}\left(\frac{1}{1+(t+\lambda)^2} +\frac{1}{1+(t-\lambda)^2}\right) \right]^{{-}2}\mathrm{d}t. \end{equation}
In the numerical results of § 4.2, the integral in (3.29) is evaluated with a symbolic programming language. Note that the Stokes lines depicted in figure 4 are not truly vertical, and are slightly curved such that they intersect the free surface at the points 
$-\lambda ^*$ and 
$\lambda ^*$. Thus, the range of integration in (3.29) should actually lie between 
$-\lambda ^*< t<\lambda ^*$; however, since 
$\lambda ^*$ is very close in value to 
$\lambda$ (for 
$\lambda =8$ and 
$\varGamma _c=0.3$, 
$\lambda ^*\approx 7.99998$), this subtlety has been ignored.
Comparisons between the analytical prediction of 
$F_k$ from (3.28) and numerical results are performed in § 4.2.
4. Numerical results
We begin in § 4.1 by verifying with numerical results our analytical predictions for the exponentially small scaling as 
$F \to 0$ for the case of a single vortex. This is given by the singulant function, 
$\chi$, from (3.9), and comparisons are made for a range of values of the vorticity, 
$\varGamma _{{c}}$. The analytical predictions of the Froude numbers for trapped waves between two point vortices, given in (3.28), are then compared to numerical predictions in § 4.2.
A detailed description of the numerical method used is given by Forbes (Reference Forbes1985), which we will briefly summarise here.
(i) The real-valued domain,
$s$, is truncated to lie between the values of $s_{{L}}$ and $s_{{R}}$. A total of $N$ discretisation points are used, such that the numerical domain is given by $s_k=s_{{L}}+(k-1)(s_{{R}}-s_{{L}})/(N-1)$ for $1 \leq k \leq N$. The unknown solution is taken to be $y^{\prime }(s)$, which we define at each grid point by $y^{\prime }_k=y^{\prime }(s_k)$. The radiation conditions are imposed by enforcing $y_1=0$, $y_1^{\prime }=0$, $x_1^{\prime }=1$, $\phi _1^{\prime }=1$, $x_1=s_{{l}}$ and $\phi _1=s_{{L}}$, and the initial guess for $y_k^{\prime }$ is either zero or a previously computed solution.(ii) Since we assume that
$y^{\prime }_k$ is known at the next grid point, the arclength relation (2.1b) yields $x^{\prime }_k$. Trapezoidal-rule integration then determines values for $x_k$ and $y_k$, which we use to find $\phi _k^{\prime }$ from Bernoulli's equation (2.1a). This process is repeated for $k=2$ to $k=N$ to find function values at every grid point.(iii) The boundary-integral equation (2.1c) is evaluated at each grid point with the known values of
$x_k$, $y_k$, $\phi ^{\prime }_k$, $x^{\prime }_k$ and $y^{\prime }$. To avoid the singularity associated with the principal-value integral $\mathcal {I}[x,y,\phi ^{\prime }]$, each unknown that is not a function of the integration variable, $t$, is instead evaluated between grid points by interpolation.(iv) This yields
$N-1$ nonlinear equations from evaluating the boundary-integral equation between each grid point, $(s_k+s_{k+1})/2$, which is closed by the $N-1$ unknowns $y^{\prime }_k$ for $k=2$ to $k=N$. Solutions are found by minimising the residual through Newton iteration. For the trapped waves studied in § 4.2, we impose an additional constraint of symmetry about $s=0$ in the real-valued solution, $y(s)$, such that the Froude number, $F$, is determined as an eigenvalue.
4.1. Waves generated by a single vortex
For the numerical results presented in this section, we have used 
$N=2000$ grid points, and a domain specified by 
$s_{{L}}=-40$ and 
$s_{{R}}=40$. In computing numerical solutions for a wide range of Froude numbers, and the values of 
$\varGamma _{{c}}=\{ 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40\}$, the exponentially small scaling as 
$F \to 0$ of the high-frequency waves present for 
$s>0$ may be measured. This is shown in the semilog plot of figure 5. We see that these lines, each of which represents solutions with a different value of 
$\varGamma _{{c}}$, are straight, and thus the amplitude of these ripples is exponentially small as 
$F \to 0$. The gradient of each of these lines is expected to closely match the exponential scaling predicted analytically, given by the singulant 
$\chi$. Along the free surface, this is given by 
$\textrm {Re}[\chi ]$ from (3.9), which takes constant values.
Figure 5. The exponentially small dependence of the wave amplitude is shown (dots) for numerical results for seven different values of 
$\varGamma _{{c}}=\{0.10,0.15,0.20,0.25,0.30,0.35,0.40\}$. Solid lines represent the analytical gradient found from the real part of 
$\chi$ in (3.9). The behaviour of this gradient for different values of the vortex strength 
$\varGamma _{{c}}$ is shown in figure 6.
In figure 6, this analytical prediction is compared to the numerical values from figure 5, and good agreement is observed. Note that there are small instabilities present in the numerical solution, which decay when the truncated domain is extended; upon which we expect the numerical results to tend towards the analytical prediction shown in figure 6.
Comparison between a numerical and asymptotic solution profile is shown in figure 7 for 
$F=0.45$ and 
$\varGamma _{{c}}=0.4$. The numerical solution is determined by the scheme detailed at the beginning of § 4, with 
$N=2000$ discretisation points in the arclength, 
$-40 \leq s \leq 40$. The asymptotic solution plots 
$x(s)=x_0(s)+F^2x_1(s)+\bar {x}(s)$ against 
$y(s)=y_0(s) + F^2 y_1(s)+\bar {y}(s)$. These early-order solutions, 
$x_0$, 
$x_1$, 
$y_0$ and 
$y_1$, are specified in (3.2a–c) and (3.3). The exponentially small components, 
$\bar {x}$ and 
$\bar {y}$, are implemented from expression (3.19). This requires knowledge of the singulant, 
$\chi$, given in (3.9), the amplitude functions 
$Y=a \mathrm {i} \varPhi /\phi _0^{\prime }$ and 
$X=-y_1^{\prime }Y$ determined from 
$\varPhi$ in (3.11), and the Stokes multiplier, 
$\mathcal {S}$, given in (3.22). A real-valued asymptotic solution is obtained through evaluating the sums 
$\bar {x}\rvert _{a=1}+\bar {x}\rvert _{a=-1}$ and 
$\bar {y}\rvert _{a=1}+\bar {y}\rvert _{a=-1}$ on the real-valued domain, 
$\sigma =s$, for 
$\textrm {Im}[\sigma ]=0$. Note that, in the determination of the constant 
$\varLambda$, its magnitude, 
$\lvert \varLambda \rvert$, has been fitted to equal that found from the corresponding numerical solution, and its argument (corresponding to a phase shift of the resultant wave) is determined from relation (3.13a,b) as 
$\textrm {arg}[\varLambda ]=a{\rm \pi} /2$.
Figure 7. For 
$F=0.45$ and 
$\varGamma _{{c}}=0.4$, a numerical solution (dashed line) is compared to an analytical solution (solid line) determined in § 3.
4.2. Trapped gravity waves between two vortices
We considered the case of two submerged point vortices analytically in § 3.5. When each vortex had the same non-dimensional circulation, 
$\varGamma _{{c}}$, and depth equal to unity, trapped waves were seen to occur for certain discrete values of the Froude number, 
$F_k$. In this section, we compare the analytical prediction for 
$F_k$ from (3.28) with numerical results. These trapped numerical solutions are found with the method detailed at the beginning of § 4. In imposing the additional constraint of symmetry to eliminate waves downstream of the vortices, the special Froude number, 
$F_k$, is determined as an eigenvalue. These results were performed for 
$N=4000$ grid points, a domain between 
$s_{{L}}=-60$ and 
$s_{{R}}=60$, and horizontal vortex placement specified as 
$\lambda =8$.
In figure 8, we plot the tail amplitude (for 
$s>\lambda$) of the asymptotic solutions for the values of 
$0.3< F<0.5$, 
$\varGamma _{{c}}=0.3$ and 
$\lambda =8$. This amplitude is equal to zero at the values of 
$F_k$ from (3.28). Figure 8 also contains additional markers denoted by (a), where 
$F = 0.3383$, and (b), where 
$F = 0.4270$. This corresponds to figure 9, where we compare numerical solutions obtained in this section, and asymptotic solutions from § 3 for those given values of 
$F$. The fit is excellent and the corresponding curves are nearly visually indistinguishable at the scale of the graphic. Finally, in figure 10, we compare the values of 
$F_k$ obtained analytically and numerically. The straight lines are the analytical prediction from (3.28), and dots represent the numerical values for 
$F_k$. The agreement between these is good for moderately small values of 
$F$; we notice that the error (horizontal distance) decreases as 
$F$ decreases from 
$0.6$ to around 
$0.4$.
Figure 8. The amplitude of oscillations present for 
$s>\lambda$ in the asymptotic solutions is shown against the Froude number, 
$F$. Here, 
$\varGamma _{{c}}=0.3$ and 
$\lambda =8$. This amplitude is equal to zero at the locations 
$F_k$ derived in (3.28). The two points marked (a) and (b) correspond to the profiles shown in figure 9.
Figure 9. Two different trapped wave solutions are shown for 
$\varGamma _{{c}}=0.3$ and 
$\lambda =8$ corresponding to (a) 
$F = 0.3383$ and (b) 
$F = 0.4270$ (as marked in figure 8). Asymptotic solutions (solid line) are compared to numerical solutions (dashed line) for (a) 
$k=22$ and (b) 
$k=14$. In each panel, the two curves are nearly indistinguishable to visual accuracy.
Figure 10. Values of the Froude number, 
$F_k$, for which the waves are trapped between each submerged vortex are shown. The numerical results of § 4.2 are shown by dots, and the analytical results from (3.28) are shown with lines. Here, 
$\lambda =8$, and for the numerical solutions 
$N=4000$, 
$s_{{L}}=-60$ and 
$s_{{R}}=60$.
Numerical solutions of the double-vortex problem are particularly challenging, however, due to the finite-difference nature of the numerical scheme, and the truncation of the infinite domain. This seems to manifest in larger (relative) errors, particularly for 
$F \lessapprox 0.4$. We have verified that increasing the number of mesh points and increasing the domain size diminishes the numerical error. In addition, we have used the linear prediction of 
$F_k=\sqrt {\lambda /(k{\rm \pi} +{\rm \pi} /2)}$ to verify our analytical results for small values of 
$\varGamma _{{c}}$. Indeed, these numerical challenges are common in exponential asymptotics comparisons, where it is necessary to use a sufficiently small value of the asymptotic parameter, but not too small so that numerical error overtakes the exponentially small character of the solution.
5. Conclusion
We have shown, through both numerical and analytical investigations, that the waves generated by submerged point vortices are exponentially small in the low-speed limit of 
$F \to 0$. Furthermore, when two submerged vortices are considered, oscillatory waves vanish downstream for certain values of the Froude number, 
$F$. Through the techniques of exponential asymptotics, we have demonstrated how these values may be derived. Their prediction relies on the understanding of singularities in the analytically continued domain that generate a divergent asymptotic expansion. The remainder to this series is exponentially small as 
$F \to 0$, and the study of the associated Stokes phenomenon yields discrete values of 
$F$ for which the waves are trapped between each vortex.
6. Discussion
The work presented here forms a basis for a number of interesting extensions involving exponentially small water waves with gravity, capillarity and/or vorticity providing singular perturbative effects.
Firstly, it should be remarked that the classical exponential asymptotics theories by, for example, Chapman & Vanden-Broeck (Reference Chapman and Vanden-Broeck2002, Reference Chapman and Vanden-Broeck2006) for capillary- and gravity-driven surface waves produced in flows over topographies rely upon the existence of closed-form conformal maps. In such problems, the governing equations for the free surface can be written in terms of a single complex-valued unknown (e.g. the complex velocity), with the velocity potential serving as the independent variable. This includes situations such as flows past polygonal boundaries (related to the availability of the Schwarz–Christoffel mapping). The arclength formulation we have used in this work provides a more general setting for wave–structure interactions with arbitrary bodies, including, for instance, flows past smoothed bodies specified in 
$(x,y)$ coordinates. Here, we have demonstrated that the exponential asymptotics can be generalised to such formulations. We expect that many of the interesting wave–structure interactions studied by Holmes et al. (Reference Holmes, Hocking, Forbes and Baillard2013) (symmetric bottom topography), Hocking et al. (Reference Hocking, Holmes and Forbes2013) (submerged semi-ellipse), and Elcrat & Miller (Reference Elcrat and Miller2006) (submerged point vortex with lower topography), for example, can be attacked using the technology we have developed here.
Secondly, the phenomenon of trapped waves is an interesting one. The exponential asymptotics interpretation, whereby waves switched on at one location (the Stokes line intersection) must be switched off at another, provides an intuitive explanation for how trapped waves form in singularly perturbative limits. The context, in our problem, relates to vortices fixed within the fluid for modelling submerged obstructions, such as the submerged cylinders studied numerically by Tuck & Scullen (Reference Tuck and Scullen1998). However, trapped waves have been detected numerically in other geometries, including submerged bumps (Hocking et al. Reference Hocking, Holmes and Forbes2013), a semi-ellipse (Holmes et al. Reference Holmes, Hocking, Forbes and Baillard2013), a trigonometric profile (Dias & Vanden-Broeck Reference Dias and Vanden-Broeck2004), spikes (Binder, Vanden-Broeck & Dias Reference Binder, Vanden-Broeck and Dias2005) and a rectangular bump (Lustri, McCue & Binder Reference Lustri, McCue and Binder2012). We expect that the ‘selection mechanism’ that produces the countably infinite set of values (3.28) is a kind of universality in eigenvalue problems (cf. Chapman et al. (Reference Chapman, Dallaston, Kalliadasis, Trinh and Witelski2022) for further discussion and examples). In our paper, we have mentioned some of the numerical challenges in verifying beyond-all-orders predictions, particularly in connection with finite-difference formulae and truncation of the infinite domain. In addition, we note that numerical verification may be further complicated by the fact that finite-difference discretisation introduces exponentially small errors in the numerical solution, as shown by Moston-Duggan, Porter & Lustri (Reference Moston-Duggan, Porter and Lustri2023).
Finally, we note that, in this paper, the forcing mechanism producing the waves was via the complex-plane singularities associated with the point vortices – then, we found that the waves were singularly perturbed due to the inertial term in Bernoulli's equation, thus producing exponentially small waves, scaling as 
$\exp (-\textrm {const.}/F^2)$. Recently, analytical solutions have been developed for pure-vorticity-driven water waves, notably in the works by Crowdy & Nelson (Reference Crowdy and Nelson2010), Crowdy & Roenby (Reference Crowdy and Roenby2014) and Crowdy (Reference Crowdy2023). In essence, we believe these solutions can serve as leading-order approximations in the regime of small surface tension; it might be expected that exponentially small parasitic ripples then exist on the surface of such vorticity-driven profiles. This would then be similar to the work of Shelton, Milewski & Trinh (Reference Shelton, Milewski and Trinh2021) and Shelton & Trinh (Reference Shelton and Trinh2022) for parasitic capillary ripples on steep gravity waves. Numerical and analytical work on this class of problems is ongoing.
Acknowledgements
We thank the anonymous referees for their helpful comments regarding the clarity of this work. We are grateful for many stimulating and motivating discussions that took place during the recent London Mathematical Society–Bath symposium ‘New Directions in Water Waves’ held at the University of Bath in July 2022. Finally, we would like to thank the Isaac Newton Institute for Mathematical Sciences, Cambridge, for support and hospitality during the programme Applicable Resurgent Asymptotics, where some work on this paper was undertaken.
Funding
J.S. and P.H.T. gratefully acknowledge support by the Engineering and Physical Sciences Research Council (EPSRC grant no. EP/V012479/1); and for work undertaken at the Isaac Newton Institute (EPSRC grant no. EP/R014604/1).
Declaration of interests
The authors report no conflict of interest.
Appendix A. Inner analysis at the singularities $\sigma =\pm \mathrm {i}$
In order to determine the constant of integration of the amplitude function 
$\varPhi (\sigma )$ from (3.10), knowledge of the inner solutions at the singularities 
$\sigma = \mathrm {i}$ and 
$\sigma =-\mathrm {i}$ is required. In this section, we study the inner boundary layer at both of these locations, for which matching with the inner limit of the outer solutions determines the constant of integration.
First, we note that, in the outer region, where 
$\sigma =O(1)$, the asymptotic series first reorder whenever
\begin{equation} \phi_0^{\prime}(\sigma) \sim F^2 \phi_1^{\prime}(\sigma), \quad y_1(\sigma) \sim F^2 y_2(\sigma), \quad x_2(\sigma) \sim F^2 x_3(\sigma). \end{equation}
In substituting for the early orders of the asymptotic solutions specified in (3.2a–c), (3.3) and (3.5), we see that each of equations (A1a–c) reorder in a boundary layer of the same width, given by 
$\sigma -a \mathrm {i} =O(F^{2/3})$. We thus introduce the inner variable, 
$\hat {\sigma }$, by the relation
\begin{equation} \sigma-a \mathrm{i}=\hat{\sigma}F^{2/3}, \end{equation}
for which 
$\hat {\sigma }=O(1)$ in the inner region. Since the asymptotic series each reorder near the two locations of 
$\sigma =\mathrm {i}$ and 
$\sigma =-\mathrm {i}$, we have again used the notation 
$a=\pm 1$ to distinguish between these two cases.
Next, to determine the form of the inner solutions, we take the inner limit of the outer series expansions for 
$\phi ^{\prime }$, 
$x$ and 
$y$, by substituting for the inner variable 
$\hat {\sigma }$ defined in (A2) and expanding as 
$F \to 0$. This yields
\begin{equation} \phi^{\prime} \sim \frac{1}{F^{2/3}} \left[-\frac{a \mathrm{i} \varGamma_{{c}}}{2{\rm \pi}} \frac{1}{\hat{\sigma}} +\cdots\right], \quad y \sim F^{2/3} \left[\frac{\varGamma_{{c}}^2}{8{\rm \pi}^2} \frac{1}{\hat{\sigma}^2}+\cdots\right], \quad x \sim a \mathrm{i} + F^{2/3} [\hat{\sigma}+\cdots], \end{equation}
where the omitted terms, represented by 
$\cdots$, are from the inner limit of lower-order terms of the outer asymptotic expansion. For instance, the next term in the inner limit of 
$\phi ^{\prime }$ is of 
$O(F^{-2/3}\hat {\sigma }^{-4})$. The form of the inner limits in (A3a–c) motivates our definition of the inner solutions, 
$\hat {\phi }(\hat {\sigma })$, 
$\hat {y}(\hat {\sigma })$ and 
$\hat {x}(\hat {\sigma })$, through the equations
\begin{equation} \phi^{\prime}={-}\frac{a \mathrm{i} \varGamma_{{c}}}{2 {\rm \pi}F^{2/3}} \frac{\hat{\phi}(\hat{\sigma})}{\hat{\sigma}}, \quad y=\frac{\varGamma_{{c}}^2F^{2/3}}{8{\rm \pi}^2} \frac{\hat{y}(\hat{\sigma})}{\hat{\sigma}^2},\quad x= a \mathrm{i} + \hat{\sigma}F^{2/3}\hat{x}(\hat{\sigma}). \end{equation}
The form of the inner variables introduced in (A4a–c) ensures that the first term in the series expansion for their outer limit will be equal to unity. Furthermore, based on the form of the inner limit of the singulant, 
$\chi$, from (3.8),
\begin{equation} \chi \sim{-}\frac{4 a \mathrm{i} {\rm \pi}^2}{3\varGamma_{{c}}^2} \hat{\sigma}^3 F^2, \end{equation}
the outer limit of the inner solutions will be a series expansion in inverse powers of 
$-4 a \mathrm {i} {\rm \pi}^2 \hat {\sigma }^3/(3 \varGamma _{{c}}^2)$. We thus introduce the variable 
$z$, defined by
\begin{equation} z={-}\frac{4 a \mathrm{i} {\rm \pi}^2}{3\varGamma^2_{{c}}} \hat{\sigma}^3 \end{equation}
to ensure that these series expansions are in inverse powers of 
$z$ alone.
A.1. Inner equation
The leading-order inner equations, as 
$F \to 0$, may now be derived by substituting (A4a–c) for the inner variables into the outer equations (2.5a)–(2.5c), yielding
\begin{gather} \hat{y}-\hat{\phi}^2=0, \end{gather}
\begin{gather}(\hat{x} + 3 z \hat{x}^{\prime})^2 -\left(\frac{1}{3z}\hat{y}-\frac{1}{2}\hat{y}^{\prime}\right)^2=1, \end{gather}
\begin{gather}\hat{\phi}\left(\hat{x} -\frac{1}{6 z}\hat{y}\right) \left(\hat{x}+3z\hat{x}^{\prime}-\frac{1}{3z}\hat{y} +\frac{1}{2}\hat{y}^{\prime} \right)=1. \end{gather}
The inner solutions, 
$\hat {\phi }(z)$, 
$\hat {y}(z)$ and 
$\hat {x}(z)$, will satisfy (A7a)–(A7c). Rather than solve these inner equations exactly, knowledge of the inner solutions is only required under the outer limit of 
$z \to \infty$ in order to match with the inner limit of the outer solutions to determine their divergent form. Thus, we will consider the following series expansions for these inner unknowns,
\begin{equation} \hat{\phi}(z)= \sum_{n=0}^{\infty} \frac{\hat{\phi}_n}{z^n}, \quad \hat{y}(z) = \sum_{n=0}^{\infty}\frac{\hat{y}_n}{z^n}, \quad \hat{x}(z) = \sum_{n=0}^{\infty}\frac{\hat{x}_n}{z^n}, \end{equation}
which hold as 
$z \to \infty$.
At leading order as 
$z \to \infty$, we have, by the definition on the inner solutions in (A4a–c),
\begin{equation} \hat{\phi}_0 =1, \quad \hat{y}_0 =1, \quad \hat{x}_0 =1. \end{equation}
Determination of 
$\hat {\phi }_n$, 
$\hat {y}_n$ and 
$\hat {x}_n$, as 
$n \to \infty$, requires the evaluation of a recurrence relation, which is now given. Firstly, substitution of expansions (A8a–c) into the inner equation (A7b) yields
\begin{equation} \left. \begin{gathered} \hat{x}_1=0, \\ \hat{x}_n = \sum_{m=1}^{n-1}\frac{(3m-1)(3n-3m-1)}{2(1-3n)} \left(\frac{\hat{y}_{m-1}\hat{y}_{n-m-1}}{36}-\hat{x}_{m}\hat{x}_{n-m} \right) \quad \text{for }n \geq 2. \end{gathered} \right\} \end{equation}Next, we substitute the same expansions into the inner equation (A7c), yielding
\begin{align} \hat{\phi}_1=\frac{1}{2}, \quad \hat{\phi}_n & =\frac{1}{36}\sum_{m=2}^{n}\sum_{q=1}^{m-1} (3q-1) (6\hat{x}_q+\hat{y}_{q-1})(6\hat{x}_{m-q}-\hat{y}_{m-q-1})\hat{\phi}_{n-m} \nonumber\\ &\quad-\frac{1}{2}\sum_{m=1}^{n} ((4-6m)\hat{x}_m -m\hat{y}_{m-1})\hat{\phi}_{n-m} \quad \text{for }n \geq 2. \end{align}
Lastly, a recurrence relation for 
$\hat {y}_n$ is found from (A7a) to be
\begin{equation} \hat{y}_1=1, \quad \hat{y}_n = \sum_{m=0}^{n}\hat{\phi}_m \hat{\phi}_{n-m} \quad \text{for }n \geq 2. \end{equation}
Assuming that 
$\hat {\phi }_{n-1}$, 
$\hat {y}_{n-1}$ and 
$\hat {x}_{n-1}$ are known, 
$\hat {x}_n$ can be determined from (A10a), which then yields a value for 
$\hat {\phi }_n$ from (A10b). Lastly, 
$\hat {y}_n$ is found by evaluating equation (A10c).
A.2. Matching and determination of the constant $\varLambda$
We now match the outer limit of the inner solution, 
$\hat {\phi }$, with the inner limit of the outer solution, 
$\phi ^{\prime }$. In writing the outer limit of the inner solution in outer variables, we have
\begin{equation} \phi^{\prime}= \frac{-a \mathrm{i} \varGamma_{{c}}}{2 {\rm \pi}}\sum_{n=0}^{\infty} \frac{F^{2n}\hat{\phi}_n}{\left(-\dfrac{4 a \mathrm{i} {\rm \pi}^2}{3 \varGamma_{{c}}^2}\right)^n(\sigma-a \mathrm{i})^{3n+1}}, \end{equation}and for the inner limit of the outer solution,
\begin{align} \phi^{\prime}&= \sum_{n=0}^{\infty} F^{2n} \phi_n^{\prime} \sim \sum_{n=0}^{\infty} -F^{2n} \chi^{\prime} \varPhi \frac{\varGamma(n+\alpha+1)}{\chi^{n+\alpha+1}}\nonumber\\ &\sim \sum_{n=0}^{\infty}-\frac{4 {\rm \pi}^2 \varLambda} {\varGamma_{{c}}^2({-}a\mathrm{i})^{1/2}}\,\mathrm{e}^{\mathcal{P}(a \mathrm{i})} \frac{F^{2n}\varGamma(n+\alpha+1)}{\left(-\dfrac{4 a \mathrm{i} {\rm \pi}^2} {3\varGamma_{{c}}^2}\right)^{n+\alpha+1}(\sigma-a \mathrm{i})^{3n+3\alpha-1/2}}. \end{align}
In the above, the inner limit of the amplitude function 
$\varPhi$ from (3.11) has been taken by defining
\begin{equation} \mathcal{P}(\sigma)=\int_{0}^{\sigma}\left[ \frac{a \mathrm{i} \phi_1^{\prime}(t)} {[\phi_0^{\prime}(t)]^3}-\frac{3}{2(t-a \mathrm{i})}\right] \,\mathrm{d}t, \end{equation}
such that 
$\mathcal {P}(\sigma )=O(1)$ as 
$\sigma \to a \mathrm {i}$. Matching (A11) with (A12), and substituting for 
$\alpha =1/2$ from (3.13a,b), determines the constant, 
$\varLambda$, as
\begin{equation} \varLambda ={-}\frac{a \mathrm{i} \mathrm{e}^{-\mathcal{P}(a \mathrm{i})}}{3 \sqrt{3}} \lim_{n \to \infty}\left(\frac{\hat{\phi}_n}{\varGamma(n+\alpha+1)}\right). \end{equation}
