2.1. The $(q,\theta )$ formulation

In this section, we repose the two-dimensional governing system (2.1) as a one-dimensional boundary-integral formulation in terms of the free-surface speed and angle. Following the traditional treatment of potential free-surface flows, we introduce the complex potential $f=\phi +\mathrm {i} \psi$. Rather than consider $f = f(z)$, we instead consider $z = z(f)$ and, hence, the flow region is known in the potential plane. The complex potential plane is shown in figure 4. From this definition, the complex velocity can be found to be ${\rm d}f/{\rm d}z = u-\mathrm {i} v$, where $(u,v)$ are the horizontal and vertical velocities.

Introducing $q$ as the streamline speed and $\theta$ as the streamline angle by the relationship $q \mathrm {e}^{-\mathrm {i} \theta }=u-\mathrm {i} v$ then yields

(2.3)\begin{equation} \frac{{{\rm d}}{} f}{{{\rm d}}{} z}=q \mathrm{e}^{-{\mathrm{i}} \theta}. \end{equation}

In this form, Bernoulli's equation (2.1c) is written as

(2.4a)\begin{equation} F^2 q^2 \frac{{{\rm d}}{} q}{{{\rm d}}{} \phi}+ \sin{(\theta)}-Bq\frac{{{\rm d}}{}}{{{\rm d}}{} \phi}\left(q\frac{{{\rm d}}{} \theta}{{{\rm d}}{} \phi}\right)=0. \end{equation}

By the analyticity of $\log {q}-\mathrm {i} \theta$, we introduce the boundary-integral equation which relates $q$ to the Hilbert transform of $\theta$ operating over the free surface. For our periodic domain from $-1/2$ to $1/2$, we integrate $\log {q} - \mathrm {i}\theta$ using Cauchy's theorem and use the periodicity conditions

(2.4b)\begin{equation} q \left(\phi -\tfrac{1}{2} \right) = q \left(\phi + \tfrac{1}{2} \right) \quad \text{and} \quad \theta \left(\phi - \tfrac{1}{2} \right) = \theta \left(\phi+\tfrac{1}{2} \right), \end{equation}

which follow from (2.1e), and the deep-water conditions (2.1d) to derive the periodic Hilbert transform given by

(2.4c)\begin{equation} \log{(q)}=\mathscr{H}[\theta](\phi)={\int\hskip -1,05em -\,}_{-{1}/{2}}^{{1}/{2}}\theta(\phi') \cot{[{\rm \pi}(\phi'-\phi)]} {{\rm d}}{} \phi'. \end{equation}

In the above, is the Cauchy principal-value integral. The above provides the crucial relationship between the components $q$ and $\theta$, and further details on the derivation of the boundary-integral relations can be found in chapter 6 of Vanden-Broeck (Reference Vanden-Broeck2010).

Finally, the energy expression (2.2) is also considered in terms of $(q,\theta )$. Noting that $x_{\phi }=q^{-1}\cos (\theta )$ and $y_{\phi }=q^{-1}\sin {\theta }$, we substitute $y=(F^2/2)(1-q^2)+Bq \theta _{\phi }$ from Bernoulli's equation to find

(2.4d)\begin{equation} \mathscr{E} = \frac{1}{E_{{hw}}}\int_{-{1}/{2}}^{{1}/{2}} \left[ \mathcal{G}_0(\phi)+B\mathcal{G}_1(\phi) +B^2\mathcal{G}_2(\phi) \right] {{\rm d}}{} \phi, \end{equation}

where we have defined components

(2.5)\begin{equation} \left.\begin{gathered} \mathcal{G}_0(\phi) = \frac{F^4}{8q}(1-q^2)(3 \cos{\theta}-2q-q^2 \cos{\theta}), \\ \mathcal{G}_1(\phi)= \frac{(1-\cos\theta)}{q}+\frac{F^2 \theta_{\phi}}{2}(2\cos{\theta}-q-q^2\cos{\theta}), \\ \mathcal{G}_2(\phi)= \frac{q \theta_{\phi}^2\cos{\theta}}{2}. \end{gathered}\right\} \end{equation}

In summary, the water-wave problem, as formulated for $q$ and $\theta$, involves the solution of (2.4a)–(2.4d). Note that the above sets of equations all involve the evaluation of $q$ and $\theta$ on the streamline $\psi = 0$.

2.2. Analytic continuation

As we show, the exponential asymptotics procedure of § 7 requires the continuation of the free-surface solutions, $q(\phi + 0\mathrm {i})$ and $\theta (\phi + 0\mathrm {i})$, into the complex plane, where $\phi \in \mathbb {C}$. This free-surface continuation procedure is depicted in figure 6. Hence, we analytically continue Bernoulli's equation (2.4a) and the boundary-integral equation (2.4c) into the complex $\phi$-plane. The independent variable $\phi$ is complexified by considering $\phi \mapsto \phi _c\in \mathbb {C}$ and, hence, $q$ and $\theta$ are analytically continued. For convenience, we relabel $\phi _c$ as $f$. Thus, Bernoulli's equation remains in an identical form to (2.4a), but with the variable $\phi$ replaced by $f$.

Figure 6.

A schematic of our analytic continuation procedure demonstrates the difference between the physical $\phi + \mathrm {i} \psi$ plane and our complexified $\phi _c$ space. The location of the principle upper- and lower-half singularities at $f^*$ and $-f^*$ of the leading-order flow field are shown by circles, and the main Stokes line from § 7 is shown dashed.

For the boundary-integral equation (2.4c), we must consider the complexification of the Hilbert transform. Let us write

(2.6)\begin{equation} \mathscr{H}[\theta]=\hat{\mathscr{H}}[\theta]-a\mathrm{i} \theta, \end{equation}

where $\hat {\mathscr {H}}[\theta ]$ is the complex-valued Hilbert transform,

(2.7)\begin{equation} \hat{\mathscr{H}}[\theta](f)=\int_{-{1}/{2}}^{{1}/{2}}\theta(\phi') \cot{[{\rm \pi}(\phi'-f)]} \, {\rm d}\phi'. \end{equation}

Note that the integral above is only evaluated along the physical free surface, parameterised in terms of $\phi ^{\prime }$, where $\theta$ takes real values.

In (2.6), we have also introduced the parameter, $a$, which is defined by

(2.8)\begin{equation} a=\left\{ \begin{array}{@{}ll} +1 & \text{for }\operatorname{Im} (f) >0,\\ -1 & \text{for }\operatorname{Im} (f) <0. \end{array} \right. \end{equation}

When the Hilbert transform relationship is extended into the upper-half-$f$-plane, $a = 1$, whereas $a = -1$ for continuation into the lower-half-$f$-plane. The validity of (2.6) as a legitimate complexification of the Hilbert transform is verified by taking $\operatorname {Im}(f) \to 0$ on the right-hand side. Then $\hat {\mathscr {H}}[\theta ]$ yields a principal value integral and residue. The residue contribution changes sign between $\operatorname {Im}(f) \to 0^+$ and $\operatorname {Im}(f) \to 0^-$, yielding the constant $a$.

In summary, the governing equations for the analytically continued $q$ and $\theta$ values are given by

(2.9a)\begin{gather} F^2 q^2 q'+ \sin{(\theta)}-Bq (q \theta' )'=0, \end{gather}

(2.9b)\begin{gather}\log(q)+a{\mathrm{i}}\theta=\hat{\mathscr{H}}[\theta], \end{gather}

(2.9c)\begin{gather}\mathscr{E} = \frac{1}{E_{{hw}}}\int_{-{1}/{2}}^{{1}/{2}} \left[ \mathcal{G}_0(\phi)+B\mathcal{G}_1(\phi) +B^2\mathcal{G}_2(\phi) \right] {{\rm d}}{} \phi, \end{gather}

(2.9d)\begin{gather}q\left(-\tfrac{1}{2}\right)=q\left(\tfrac{1}{2}\right) \quad \text{and} \quad q^{\prime}\left(-\tfrac{1}{2}\right)=q^{\prime}\left(\tfrac{1}{2}\right). \end{gather}

Note that although (2.9a) and (2.9b) are evaluated through complex $f$-space, the energy condition is most easily evaluated on the physical free surface. Here and henceforth, we use primes ($\prime$) to denote differentiation in $f$. This system will be solved in § 4 with an expansion holding under the limit of $B \to 0$.

3.1. Asymptotic inconsistencies in Longuet-Higgins (1963)

Numerical evidence was provided by Shelton et al. (Reference Shelton, Milewski and Trinh2021) (see figure 3) to demonstrate that, for those solutions exhibiting small-scale ripples on an underlying gravity wave, the amplitude of these parasitic ripples is exponentially small as $T' \to 0$. Solutions that display such exponentially small behaviour cannot be described purely by a typical Poincaré expansion which contains only algebraic powers of the small parameter; their description will instead appear beyond all orders of the standard Poincaré expansion.

We now review Longuet-Higgins’ approach in our non-dimensional formulation (using the Bond number, $B$, and Froude number, $F$, in (1.2a,b) instead of $T'$ and $c$). We start with the integrated form of Bernoulli's equation from (1.1) given in terms of $y$ and the streamline speed, $q$, as

(3.3)\begin{equation} \frac{F^2}{2}q^2+y-B\frac{\partial q}{\partial \psi}=\text{const}., \end{equation}

where the derivative in the $\psi$ direction can be converted into a derivative the $\phi$ direction via the Cauchy–Riemann equations. In his §3, Longuet-Higgins considered a perturbation $(\bar {y},\bar {q})$ about the gravity-wave $(y_0,q_0)$ with the truncations from (3.1a–c) to find

(3.4)\begin{equation} \frac{F^2}{2}(q_0^2+2q_0\bar{q}+\bar{q}^2)+(y_0+\bar{y})-B\left(\frac{\partial q_0}{\partial \psi} +\frac{\partial \bar{q}}{\partial \psi}\right)=\text{const}. \end{equation}

Here, the $O(1)$ terms, $F^2 q_0^2/2+ y_{0}=\textrm {const.}$, are satisfied exactly as this is the gravity-wave equation with solutions $(y_0,q_0)$. Thus, we obtain

(3.5)\begin{equation} \underbrace{F^2q_0\bar{q}+\bar{y}-B \frac{\partial q_0}{\partial \psi}}_{O(B)} -\underbrace{B\frac{\partial \bar{q}}{\partial \psi}}_{O(B^2)}={-}\underbrace{\frac{F^2\bar{q}^2}{2}}_{O(B^2)}. \end{equation}

The asymptotic behaviour indicated by the under-braced quantities follows by making the standard assumption that the leading corrections, $\bar {y}$ and $\bar {q}$, are both of $O(B)$. Consequently, $\bar {q} \ll q_0$, and so Longuet-Higgins neglected the nonlinear term $\bar {q}^2$ on the right-hand side of this equation. However, the $O(B^2)$ term on the left-hand side was not neglected. This assumption, which appears in his (5.1), is asymptotically inconsistent. In fact, this inconsistency is how Longuet-Higgins was able to produce approximations to an a priori exponentially small capillary ripple, because, otherwise, all corrections are ripple-free and algebraic in $B$.

The above asymptotic inconsistency is somewhat typical in early models of many exponential asymptotic problems. There are two (formally correct) methods to proceed with (3.5).

1. (i) We may correctly treat $\bar {y}$ and $\bar {q}$ to both be of $O(B)$. The leading-order terms in (3.5) are, thus,

   (3.6)\begin{equation} F^2q_0\bar{q}+\bar{y}-B \frac{\partial q_0}{\partial \psi}=0, \end{equation}
   and would yield the $O(B)$ capillary correction term. The procedure could be continued to quadratic orders of $B$ and higher, but the resultant perturbative solution would never yield an exponentially small ripple. In essence, this is a derivation of the regular perturbative expansion and leads to the analysis of § 4.

2. (ii) Alternatively, we may consider $\bar {y}$ and $\bar {q}$ to both scale as $\sim \mathrm {e}^{-\alpha /B}$, i.e. for solutions to be of Wentzel–Kramers–Brillouin–Jeffreys (WKBJ) type. As differentiation of this ansatz yields a factor of $1/B$, the dominant terms in (3.5) change to

   (3.7)\begin{equation} \underbrace{F^2 q_0 \bar{q} - B \frac{\partial \bar{q}}{\partial \psi}}_{O(\mathrm{e}^{-\alpha/B})} = \underbrace{B q_0'}_{O(B)}. \end{equation}
   The form of the above equation would allow for the correct prediction of the WKBJ phase, $\alpha$, but not the correct prefactor (amplitude); this is on account of the fact the right-hand side is the result of a one-term truncation of the Poincaré expansion (3.1a–c). Instead, the correct procedure must involve additional terms of the regular expansion. In general, the right-hand side is of $O(B^N)$ with $N \to \infty$ as $B \to 0$. In order to derive the exponentially small ripples, we must optimally truncate with $N$ chosen carefully (Chapman et al. Reference Chapman, King and Adams1998).

Longuet-Higgins had worked with the asymptotically inconsistent (3.5), with the right-hand side set to zero, and this was used to derive the solution (3.2).

As shown in § 7.3, the ripples have the analytical behaviour

(3.8)\begin{equation} q_{{exp}}(\phi)=\varLambda \mathcal{F}(\phi)\exp{\left(-\frac{{\chi}(\phi)}{B}\right)}, \end{equation}

where $\varLambda$ is a constant coefficient, $\mathcal {F}(\phi )$ is a functional prefactor and $\chi (\phi )$ is the exponentially small dependence of the solution, which is related to the quantity $\alpha (\phi )$. These components will be significantly different than those derived by Longuet-Higgins in (3.2). In order to be correct, the above expression must be derived through optimal truncation of the standard asymptotic expansion, rather than using the one-term truncation in (3.1a–c).

We note that it is still nevertheless possible to capture exponentially small behaviour with the truncation (3.1a–c) used by Longuet-Higgins. A comprehensive review of truncations of this type, for the case of free-surface flows, was given by Trinh (Reference Trinh2017) who, aided by the use of exponential asymptotics, discussed how the functional form of the exponentially small waves changes when different truncations are made. The type utilised here by Longuet-Higgins in (3.1a–c) is an $N=1$ truncation as only one term of the asymptotic series is included. Although this truncation (if dealt with in an asymptotically consistent manner) can predict the correct exponentially small scaling of the solution, the functional form of the prefactor and its magnitude [cf. (3.2b)] will be incorrect.

3.2. The choice of integration in the exponential argument

We now discuss the second issue with Longuet-Higgins’ analytical solution, which is that (3.2) predicts an $O(1)$ solution magnitude. For real values of $\phi$, $\alpha$ takes purely real values. Thus, as his solution contains $\mathrm {e}^{-\mathrm {i} c \alpha /T'}$, only a rapidly oscillating waveform of wavelength $O(\epsilon )$ is predicted. The issue is not precisely related to the functional form of the exponential argument, because modulo the scalings, it can be confirmed via our work that

(3.9)\begin{equation} -\frac{1}{B} \frac{{{\rm d}}{}\chi}{{{\rm d}}{}\phi} \propto{-}\frac{\mathrm{i} c}{T'} \frac{{{\rm d}}{}\alpha}{{{\rm d}}{}\phi}. \end{equation}

However, Longuet-Higgins restricts $\phi$ to take real values and forces the starting point of integration in $\alpha (\phi )$ to be at $\phi =0$. This is later matched to an ad hoc simplification near the crest of the wave. This misses a fundamental step in the determination of the parasitic ripples because, as we show, their existence is intimately connected with the singularities of ${\chi }^{\prime }{(\phi )}$ in the analytic continuation of the free surface. In order to correctly resolve the Stokes phenomenon in § 7, integration in our expression for $\chi$ must begin from such singularities, and results in a path of integration through the complex-valued domain. The final result produces a complex-valued $q_{exp}$, which is paired with a conjugate contribution to in order to produce a real-valued solution with both exponentially small phase and amplitude.

5.1. Singularities in the analytic continuation of $q_0$ and $\theta _0$

The singular scaling of $q_0$ and $\theta _0$ is now considered. We let $f^*$ denote the ‘crest’ singularity in the upper-half-$f$-plane. We leave the constant, $a$, unspecified and take the limit of $f \rightarrow af^*$. First, it can be verified a posteriori that as $f \to af^*$, $|\operatorname {Im}\theta _0| \to \infty$ and

(5.2)\begin{equation} \sin{\theta_0}=\frac{1}{2 \mathrm{i}}\left[\mathrm{e}^{\mathrm{i} \theta_0}-\mathrm{e}^{-\mathrm{i} \theta_0}\right] \sim \frac{a}{2\mathrm{i} }\mathrm{e}^{a\mathrm{i} \theta_0}. \end{equation}

We multiply Bernoulli's equation (4.2a) by $q_0$, and use the above scaling for $\sin {\theta _0}$ to find

(5.3)\begin{equation} F_0^2 q_0^3 \frac{{{\rm d}}{} q_0}{{{\rm d}}{} f}={-}q_0\sin{\theta_0} \sim{-}\frac{a}{2\mathrm{i} }q_0 \mathrm{e}^{a\mathrm{i} \theta_0}. \end{equation}

However, in taking the exponential of the boundary-integral equation (4.2b), we have

(5.4)\begin{equation} q_0 \mathrm{e}^{a \mathrm{i}\theta_0} = \mathrm{e}^{\hat{\mathscr{H}}[\theta_0]}. \end{equation}

Note that the complex Hilbert transform is applied to $\theta _0$ and integrated over the free surface, where $\theta _0 = O(1)$. Thus $q_0 \mathrm {e}^{a \mathrm {i}\theta _0}$ is also of order unity and we conclude from (5.3) that $q^3_0q_0^{\prime }$ tends to a constant as $f \to af^*$. Integration then yields the following singular behaviour for $q_0$,

(5.5)\begin{equation} q_0 \sim c_a(f-af^*)^{{1}/{4}}. \end{equation}

In addition, the scaling for $\mathrm {e}^{a\mathrm {i} \theta _0}$ is found from (5.3), giving

(5.6)\begin{equation} \mathrm{e}^{a\mathrm{i} \theta_0} \sim \frac{-a\mathrm{i} F_0^2c_a^3}{2} (f-af^*)^{-{1}/{4}}. \end{equation}

Combining these results for $q_0$ in (5.5) and $\theta _0$ in (5.6) gives the scaling for the complex velocity,

(5.7)\begin{equation} \frac{{{\rm d}}{} f}{{{\rm d}}{} z} \sim c_a\left(\frac{-a \mathrm{i} F_0^2c_a^3}{2}\right)^{{-}a}(f-af^*)^{({a+1})/{4}}. \end{equation}

Note that $a=1$ recovers the same singular behaviour of Grant (Reference Grant1973) in the upper-half plane, shown in (5.1).

5.2. The apparent paradox of a singularity in the lower-half plane

We see from (5.5) for $q_0$ that a singularity exists ‘within the fluid’ in the lower-half plane at $f=-f^{*}$. This is in contrast to the regular behaviour near the same location provided by Grant's result. Our apparent prediction of singular behaviour in the flow field can readily be resolved by noting that this singularity is for the analytically continued variable, originally relabelled from $q_c \to q$ in § 2.2. It is, thus, important to distinguish between the complexified and ‘physical’ streamline speeds $q_c$ and $q_{{phys}}$, and angles $\theta _c$ and $\theta _{{phys}}$. These physical variables are found by taking the magnitude and argument of the complex velocity $q_0 \mathrm {e}^{-\mathrm {i} \theta _0}$ as in (5.7), which is regular for $a=-1$, yielding

(5.8a,b)\begin{equation} q_{{phys}}=\left\lvert q_c \mathrm{e}^{- \mathrm{i} \theta_c}\right\rvert \quad \text{and} \quad \theta_{{phys}}= \operatorname{Arg}{ \left(q_c \mathrm{e}^{- \mathrm{i} \theta_c} \right)}. \end{equation}

Thus, as $f \to -f^*$ these physical values are regular for the leading-order Stokes wave solution. Only by recombining $q_0 \mathrm {e}^{-\mathrm {i} \theta _0}$ to find the physical values within the fluid have these singular terms cancelled out.

6.1. Determination of $\gamma$ and $\varLambda$

At this point, we have determined the key components, $Q$, $\varTheta$ and $\chi$, that appear in the factorial-over-power ansatz (6.1a,b). This leaves the value of the constants $\gamma$ and $\varLambda _a$. Note that our asymptotic series (4.1a,b) reorders as $f \to af^*$ (for which $q_0 = O(Bq_1)$ for instance) and the matched asymptotics procedure that results in investigating this limit yields $\gamma$ and $\varLambda _a$.

In order to determine the constant $\gamma$, we take the limit $f \to af^*$ and match the order of the singularity of the divergent ansatz, valid for $n$ large, to the low-order behaviour. Setting $n=0$ in (6.9) and taking the limit of $f \to af^{*}$ yields

(6.10)\begin{equation} \left.q_n \right\vert_{n = 0} = O \left(\frac{q_0^2}{\chi^{\gamma}} \exp \left(a{\mathrm{i}}\theta_0+a \mathrm{i}\int_{0}^{f}\left[ \frac{\cos{\theta_0}}{F_0^2 q_0^{3}} - F_0^2q_1-2F_0F_1q_0\right]\, {\rm d}f'\right)\right). \end{equation}

From the scalings of $q_0$ and $\theta _0$ in § 5.1, and the scaling of $q_1$ in Appendix A we find that

(6.11a) \begin{gather} \left.\begin{gathered} \chi^{\gamma} = O ((f-af^*)^{{5\gamma}/{4}}), \\ q_0^2 \exp \left( a \mathrm{i}\int_{0}^{f} \left[\frac{\cos{\theta_0}}{F_0^2 q_0^{3}} -F_0^2q_1-2F_0F_1q_0 \right]{\rm d}f \right) = O\left((f-af^*)^{{5}/{4}}\right). \end{gathered}\right\}\end{gather}

We substitute the above into (6.10) and match to $q_0 = O(f-af^*)^{1/4}$ to find

(6.12)\begin{equation} \gamma=\tfrac{4}{5}. \end{equation}

As is the case in many exponential asymptotic analyses, the determination of the constant prefactor, $\varLambda _a$, is often the most troublesome aspect of the procedure. For our purposes, it will be sufficient to know that $\varLambda _a$ is a non-zero constant, and can be determined via the solution of a numerical recursion relation. Specifically, it is found by matching the ‘inner’ limit of $q_n$ from the divergent form (6.9) with the ‘outer’ limit of the inner solution for $q$ near $f = af^*$. This analysis is performed in Appendix B, yielding

(6.13)\begin{equation} \varLambda_a={-}\frac{2 \mathrm{i} f^*}{F_0^2 c_a^4} \mathrm{e}^{-\mathcal{P}(af^*)}\left( \frac{4 a \mathrm{i} F_0^2 c_a}{5} \right)^{{4}/{5}} \lim_{n \to \infty} \frac{\hat{q}_n}{\varGamma(n+\gamma)}. \end{equation}

Here, $\hat {q}_n$ is the $n$th term of the outer limit of an inner solution holding near $f=af^*$, and can be determined by recurrence relation (B14). The constant $c_a$ is the prefactor of the singular scaling of $q_0$ from (5.6) whereas $\mathcal {P}(af^*)$ is given in (A7). We do not need to work with the precise value of $\varLambda _a$; however, later in §§ 6.2 and 7.3, the fact that $\varLambda _{1}$ and $\varLambda _{-1}$ are complex conjugates will be crucial to obtain a real-valued solution on the free surface. As the prefactor, $\varLambda _a$, only has a scaling effect on the solutions (and is independent of $B$), it will be convenient to choose a specific value for visualisation purposes in § 8.

6.2. The divergence along the free surface

In order to capture the divergence of $q_n$ along the free surface, $\textrm {Im}[f]=0$, we must include the effects of the two symmetrically placed crest singularities indexed by $a = \pm 1$. We thus write

(6.14)\begin{equation} \mathfrak{q}_{n}=q_n \rvert_{a=1}+q_n\rvert_{a={-}1}. \end{equation}

By the results of § B.3, the constants $\varLambda _{1}$ and $\varLambda _{-1}$ are the complex conjugates of one another. With regards to the two singulants, $\chi _1$ and $\chi _{-1}$, we may split the path of integration via

(6.15)\begin{equation} \chi_a(\phi)= a \mathrm{i} F_0^2 \left[\int_{a f^{*}}^{0} + \int_{0}^{\phi}\right] q_0(f') \, {{\rm d}}{} f', \end{equation}

for $f = \phi$ along the real axis. As $q_0$ takes real values on the free-surface, $\textrm {Im}[f]=0$, the second integral above is seen to take purely imaginary values. By the Schwarz reflection principle, $q_0$ evaluated on the imaginary axis between $-af^*$ and $af^*$ is purely real and symmetric about the origin. Therefore the first integral on the right-hand side of (6.15) is purely real and takes the same value regardless of the choice of $a$. Thus, $\chi _{-1}$ and $\chi _{1}$ are also the complex conjugates of one another on the free surface.

Due to this behaviour of $\varLambda _a$ and $\chi _a$, we write

(6.16a,b)\begin{equation} \varLambda_a=\lvert \varLambda_1 \rvert \exp [{a \mathrm{i} \arg{\varLambda_1}}] \quad \text{and} \quad \chi_a(\phi)= \lvert \chi_1(\phi) \rvert \exp [{a \mathrm{i} \arg{\chi_1(\phi)}}], \end{equation}

which upon substitution into $\mathfrak {q}_{n}=q_n \rvert _{a=1}+q_n\rvert _{a=-1}$ yields

(6.17)\begin{equation} \mathfrak{q}_n(\phi)= \frac{2 \lvert \varLambda_1 \rvert q_0^2 \varGamma{(n+\gamma)}}{\lvert \chi_1(\phi)\rvert^{n+\gamma}} \cos{\left[\arg{\varLambda_1}-(n+\gamma)\arg{\chi_1}(\phi)+\theta_0+I(\phi) \right]}, \end{equation}

where we have defined

(6.18)\begin{equation} I(\phi)=\int_{0}^{\phi} \left( \frac{\cos{\theta_0}}{F_0^2 q_0^{3}} -F_0^2q_1-2F_0F_1q_0 \right) {{\rm d}}{}{\phi'}. \end{equation}

Thus, the above form (6.17) captures the real-valued divergence on the free-surface.

We have successfully derived an expression for the late term divergence on the axis in (6.17) and off the axis in (6.9).

7.1. Analysis of the Stokes lines

To find the Stokes lines for our problem, we apply conditions (7.1a,b) to our expression for the singulant, $\chi _a$, given in (6.5) as

(7.2)\begin{equation} \chi_a(f) =a{\mathrm{i}}F_0^2\int_{af^*}^{f} q_0(f') \,{\rm d}f'. \end{equation}

Here, integration begins at the principal singularity, $f'=af^{*}$, that lies in the analytic continuation of the free surface. Note that unlike many traditional studies in exponential asymptotics, the determination of the singulant function requires the leading-order solution, $q_0$, for which there does not exist a closed-form analytical solution. We use numerical values of $q_0$ to evaluate the singulant, $\chi _a$.

The procedure is as follows. Given a fixed value of the energy, we obtain numerical values of $q_0$ and $\theta _0$ along the free-surface $\textrm {Im}[f]=0$ using the numerical computations of Shelton et al. (Reference Shelton, Milewski and Trinh2021) or any standard procedure for calculating gravity Stokes waves (cf. Vanden-Broeck Reference Vanden-Broeck1986). Next, the analytic continuation method of Crew & Trinh (Reference Crew and Trinh2016) is used to find $q_0$ and $\theta _0$ in the complex $f$-plane. Values for $\chi _a$ are then found across the domain by integrating $q_0$ along paths originating at either singularity. Graphs of the critical contours of $\textrm {Im}[\chi _a]$ and $\textrm {Re}[\chi _a]$ are given in figure 7 for the two choices of $a=1$ and $a=-1$. We see that there are two Stokes lines along the imaginary axis from $f=-f^*$ to $f=f^*$, one for $a=1$ and another for $a=-1$, which intersect with the free surface at the wave crest $\phi =0$.

Figure 7.

Values of the singulant, $\chi _a$, are shown from the upper-half singularity in (a), with $a = 1$, and from the lower-half singularity in (b), with $a = -1$. The Stokes lines, which satisfy conditions (7.1a,b), are shown by the thick lines in the grey-shaded regions. This configuration corresponds to the energy $\mathscr {E} = 0.3804$ for which the upper-half plane singularity is at $f^* \approx 0.07776\mathrm {i}$. The chosen branch cuts for each of these singularities are shown by dashed lines.

Note that only the Stokes lines that intersect with the free-surface, $\textrm {Im}[f]=0$, are considered; other Stokes lines would indicate a switching on or switching off of exponentials in the general complex plane, but are not associated with the physical production of surface ripples.

7.2. Optimal truncation

In order to capture the exponentially small components of the solution, which do not appear in the Poincaré series (4.1a,b), we truncate the series at $n=N-1$ by considering

(7.3a–c)\begin{equation} q= \underbrace{\sum_{n=0}^{N-1}B^n q_n}_{q_r} + \bar{q}, \quad \theta = \underbrace{\sum_{n=0}^{N-1}B^n \theta_n}_{\theta_r} + \bar{\theta} \quad \text{and} \quad F = \underbrace{\sum_{n=0}^{N-1} B^n F_n}_{F_r} + \bar{F}, \end{equation}

and, thus, we have introduced the notation of $q_r$, $\theta _r$ and $F_r$ for the truncated regular expansions of the solutions and eigenvalue.

We demonstrate that when truncated optimally at the point where two consecutive terms are of the same order, that is, choose $N$ such that $\lvert B^N q_N \rvert \sim \lvert B^{N+1}q_{N+1}\rvert$, the remainders $\bar {q}$, $\bar {\theta }$ and $\bar {F}$ will be exponentially small. This point of optimal truncation is given by

(7.4)\begin{equation} N = \frac{\lvert\chi_a\rvert}{B} + \rho, \end{equation}

where $\rho \in [0, 1)$ is a bounded number to ensure $N$ is an integer.

Substituting these into the boundary-integral equation (2.9b) yields a relationship between $\bar {\theta }$ and $\bar {q}$, given by

(7.5)\begin{equation} \bar{\theta}=\frac{a \mathrm{i}\bar{q}}{q_r}-a \mathrm{i} \xi_{{int}}-a \mathrm{i} \hat{\mathscr{H}}[\bar{\theta}]+O(\bar{q}^2). \end{equation}

Similarly we can insert the truncations (7.3a–c) into Bernoulli's equation (2.9a). This gives a second-order differential equation for $\bar {q}$ and $\bar {\theta }$. Upon substituting for $\bar {\theta }$ from (7.5), this is reduced down to an equation for $\bar {q}$ only. Furthermore, we neglect the Hilbert transform of the remainder, $\hat {\mathscr {H}}[\bar {\theta }]$, as this is anticipated to be exponentially subdominant. This yields

(7.6)\begin{align} &\left[ a \mathrm{i} Bq_r \right] \bar{q}^{\prime \prime} + \left[{-}F_r^2q_r^2-a \mathrm{i} Bq_r^{\prime}+Bq_r\theta_r^{\prime}-a \mathrm{i} Bq_r \xi_{{{int}}}^{\prime}\right] \bar{q}^{\prime} \nonumber\\ &\quad + \left[ -\frac{a \mathrm{i} \cos{\theta_r}}{q_r}-2F_r^2q_rq_r^{\prime}+\frac{a \mathrm{i} B(q_r^{\prime})^2}{q_r} +Bq_r^{\prime}\theta_r^{\prime} \right.\nonumber\\ &\quad \left.- \vphantom{\frac{a \mathrm{i} \cos{\theta_r}}{q_r}} a \mathrm{i} Bq_r^{\prime \prime}+2Bq_r\theta_r^{\prime \prime}- a \mathrm{i} B q_r^{\prime}\xi_{{int}}^{\prime}-2a \mathrm{i} Bq_r \xi_{{int}}^{\prime \prime} \right] \bar{q} -2F_rq_r^2q_r^{\prime}\bar{F} = \mathcal{R} + O(\bar{q}^2). \end{align}

This is a second-order differential equation for $\bar {q}$, in which the forcing terms on the right-hand side are of $O(B^N)$. A similar equation was derived by Trinh (Reference Trinh2017) for the low-Froude limit of gravity waves. Here, we have introduced the forcing terms $\xi _{{int}}$ and $\xi _{{bern}}$ arising from the Poincaré expansion in the boundary-integral and Bernoulli's equations as

(7.7a)\begin{gather} \xi_{{int}}=\hat{\mathscr{H}}[\theta_r]-a \mathrm{i} \theta_r -\log{q_r}, \end{gather}

(7.7b)\begin{gather}\xi_{{bern}}=F_r^2q_r^2q_r^{\prime}+\sin{(\theta_r)}-B(q_r^2\theta_r^{\prime \prime}+\theta_r^{\prime}q_r^{\prime}q_r), \end{gather}

(7.7c)\begin{gather}\mathcal{R} = \xi_{{bern}}-a \mathrm{i} \cos{\theta_r}\xi_{{int}}+a \mathrm{i} Bq_rq_r^{\prime}\xi_{{int}}^{\prime}+a \mathrm{i} Bq_r^2\xi_{{int}}^{\prime \prime}. \end{gather}

Due to the truncation at $n=N-1$, equation (7.6) is satisfied exactly for every order up to and including $B^{N-1}$ because $\xi _{{int}} = O(B^N)$ and $\xi _{{bern}} = O(B^N)$.

7.3. Stokes line smoothing

We now seek a closed-form asymptotic expression for $\bar {q}$ and the terms switched-on across Stokes lines. We start with the homogeneous form of (7.6), in which the terms on the right-hand side and $\bar {F}$ are neglected. Following the exponential asymptotics methodology established in, e.g., § 4 of Chapman & Vanden-Broeck (Reference Chapman and Vanden-Broeck2006), we note that the homogeneous problem has solutions of the form,

(7.8)\begin{equation} \bar{q}_{homog} \sim Q_a(f) \exp\left(-\frac{\chi_a(f)}{B}\right), \end{equation}

where $\chi _a(f)$ and $Q_a(f)$ satisfy those same equations as found for the late-term ansatz via (6.4) and (6.6). To observe the Stokes phenomenon and the switching of exponentials, we now include the forcing terms on the right-hand side of (7.6) for $\bar {q}$. We consider a solution of the form

(7.9)\begin{equation} \bar{q}(f)=A_a(f)Q_a(f)\exp\left(-\frac{\chi_a(f)}{B}\right), \end{equation}

where the Stokes multiplier $A_a(f)$ is introduced to capture the switching behaviour that occurs across the Stokes lines. When the truncation point, $N$, is chosen optimally as in (7.4), $\bar {q}$ will be seen to be exponentially small and will change in magnitude across the lines where $\operatorname {Im}[\chi _a]=0$ and $\operatorname {Re}[\chi _a]\geqslant 0$.

The algebra for this procedure follows very similarly to, e.g., Chapman et al. (Reference Chapman, King and Adams1998), Chapman & Vanden-Broeck (Reference Chapman and Vanden-Broeck2006) and Trinh (Reference Trinh2017). Thus, when the exponential form of (7.9) for $\bar {q}$ is substituted into (7.6), the dominant balance at leading order is identically satisfied by our choice of $\chi$ determined in (6.4). The first non-trivial balance occurs at $O(\mathrm {e}^{-\chi /B})$, which also involves the forcing terms on the right-hand side. We extract the $O(B^N)$ terms from $\mathcal {R}$ in (7.7c), and this yields $\mathcal {R} \sim -q_0^2\theta _{N-1}^{\prime \prime }B^N$. The governing equation for $A_a$ is then given by

(7.10)\begin{equation} \left[F_0^2q_0^2 Q_a \mathrm{e}^{-{\chi_a}/{B}} \right] \frac{{\rm d}A_a}{{\rm d}f} \sim{-}a \mathrm{i} q_0 q_{N-1}^{\prime \prime}B^N, \end{equation}

where we have used $\theta _{N-1}^{\prime \prime } \sim a \mathrm {i} q_0^{-1} q_{N-1}^{\prime \prime }$ from the boundary-integral equation (6.3b).

By substituting in the factorial-over-power form for $q_{N-1}^{\prime \prime }$ from (6.1a,b), and using the chain rule to change differentiation to be in terms of $\chi _a$, we find

(7.11)\begin{equation} \frac{{{\rm d}}{}{A_a}}{{{\rm d}}{}{\chi_a}}=\frac{B^N \mathrm{e}^{\chi_a/B}\varGamma(N+1+\gamma)}{\chi_a^{N+1+\gamma}}. \end{equation}

This is now of an equivalent form to that found by Chapman & Vanden-Broeck (Reference Chapman and Vanden-Broeck2006) for the low-Froude limit of gravity waves (cf. their (4.4)). In brief, the procedure is as follows. First, we write $\chi _a = r_a \mathrm {e}^{\mathrm {i} \vartheta _a}$ and truncate optimally via (7.4) with $N = r_a/B + \rho$. Examination of the differential equation (7.11) shows that there exists a boundary layer at $\vartheta _a = 0$ and indeed this is the anticipated Stokes line where $\operatorname {Im}[\chi _a]=0$. The appropriate inner variable near the Stokes line is $\vartheta _a=B^{1/2}\overline{\vartheta}_{a}$ and (7.11) can then be integrated to show

(7.12)\begin{equation} A_a(f)=C_a + \frac{\sqrt{2 {\rm \pi}} \mathrm{i}}{B^{\gamma}} \int_{-\infty}^{\overline{\vartheta}_a \sqrt{r_a}} \exp{({-}t^2 / 2)} \, {{\rm d}}{} t, \end{equation}

where $C_a$ is constant. Taking the outer limit of $\bar {\vartheta } \to \infty$, we then see that across the Stokes line, there is a jump of

(7.13)\begin{equation} A_a(\vartheta_a \to 0-)-A_a(\vartheta_a \to 0+)=\frac{2{\rm \pi} \mathrm{i} }{B^{\gamma}}. \end{equation}

As concerns the relationship between Stokes-line contributions from $f = f^*$ and $f = -f^*$, note that as $\chi _{1}$ is the complex conjugate of $\chi _{-1}$, we have $\vartheta _{1}=-\vartheta _{-1}$. Thus, we anticipate that $C_1$ switches to $C_1 + 2{\rm \pi} \mathrm {i}/B^\gamma$ as one proceeds from left-to-right across the Stokes line from $f = f^*$. This is shown in figure 8(a). On the other hand, $C_{-1}$ switches to $C_{-1} + 2{\rm \pi} \mathrm {i}/B^\gamma$ proceeding from right-to-left across the Stokes line from $f = -f^*$. This is shown in figure 8(b). We emphasise that the above Stokes smoothing procedure only provides the local change of the prefactor, $A_a$, across the Stokes line. Determination of the constant, $C_a$, will follow from imposition of the boundary conditions.

Figure 8.

The Stokes smoothing procedure is visualised in the f-plane for (a) $a=1$ and (b) $a=-1$.

Returning now to (7.9), we write the leading-order exponentials on the axis, $\textrm {Im}[f]=0$, via $\bar {\mathfrak {q}}=\bar {q}\rvert _{a=1}+\bar {q}\rvert _{a=-1}$, either as an inner solution

(7.14a)\begin{equation} \bar{\mathfrak{q}}(\phi)=A_1(\phi)Q_1(\phi)\exp\left(-\frac{\chi_1(\phi)}{B}\right)+A_{{-}1}(\phi)Q_{{-}1}(\phi)\exp\left(-\frac{\chi_{{-}1}(\phi)}{B}\right), \end{equation}

for which $A(\phi )$ is given by (7.12), or as an outer solution by

(7.14b)\begin{equation} \bar{\mathfrak{q}}(\phi) \sim \left\{\begin{array}{@{}ll} C_1 \left.\left(Q_{a} \mathrm{e}^{-{\chi_{a}}/{B}}\right)\right\vert_{a=1} + \left\{ C_{{-}1}+\dfrac{2 {\rm \pi}\mathrm{i}}{B^\gamma}\right\} \left.\left(Q_{a} \mathrm{e}^{-{\chi_{a}}/{B}}\right)\right\vert_{a={-}1} & \text{for $\phi<0$},\\ \left\{ C_1 + \dfrac{2 {\rm \pi}\mathrm{i}}{B^\gamma}\right\} \left.\left(Q_{a} \mathrm{e}^{-{\chi_{a}}/{B}}\right) \right\vert_{a=1} +C_{{-}1} \left.\left(Q_{a} \mathrm{e}^{-{\chi_{a}}/{B}}\right) \right\vert_{a=-1} & \text{for $\phi>0$}. \end{array}\right. \end{equation}

In (7.14b), the constants, $C_1$ and $C_{-1}$, will be determined by enforcing periodicity on $\bar {\mathfrak {q}}$ and $\bar {\mathfrak {q}}^{\prime }$, as given by

(7.15a,b)\begin{equation} \bar{\mathfrak{q}}\left({-}\tfrac{1}{2}\right) =\bar{\mathfrak{q}} \left(\tfrac{1}{2}\right) \quad \text{and} \quad \bar{\mathfrak{q}}^{\prime} \left({-}\tfrac{1}{2}\right) =\bar{\mathfrak{q}}^{\prime}\left(\tfrac{1}{2}\right). \end{equation}

The second relation above arose by evaluating the derivative of periodicity condition (2.4b) at $\phi =0$. In writing $C_1=C_{1}^R + \mathrm {i} C_{1}^I$ and $C_{-1}=C_{-1}^R+ \mathrm {i} C_{-1}^I$, we have four unknowns balancing the four equations from the real and imaginary parts of (7.15a,b). Using $\varLambda _a=\lvert \varLambda _1 \rvert \mathrm {e}^{a \mathrm {i} \arg {\varLambda _1}}$ from (6.16a,b) and $\chi _a=\operatorname {Re}[\chi _1]+ a \mathrm {i} \operatorname {Im}[\chi _{1}]$ then yields the solutions

(7.16)\begin{equation} \left.\begin{gathered} C_{1}^I={-}\frac{\rm \pi}{B^{\gamma}}, \quad C_{1}^R={-}\frac{\rm \pi}{B^{\gamma}}\frac{\cos{[G(1/2)]}}{\sin{[G(1/2)]}}, \\ C_{{-}1}^I={-}\frac{\rm \pi}{B^{\gamma}}, \quad C_{{-}1}^R={-}\frac{\rm \pi}{B^{\gamma}}\frac{\cos{[G(1/2])}}{\sin{[G(1/2)]}}, \end{gathered}\right\} \end{equation}

where

(7.17)\begin{equation} G(\phi)=\theta_0(\phi)+\int_{0}^{\phi}\left[ \frac{\cos{\theta_0}}{F_0^2 q_0^3}-F_0^2q_1-2F_0F_1q_0-\frac{F_0^2 q_0}{B}\right] {{\rm d}}{}{\phi}. \end{equation}

Solutions are not possible when $\sin [G(1/2)]=0$, from which we obtain the following discrete set of values of $B$,

(7.18)\begin{equation} B_n = \frac{F_0^2 \displaystyle\int_{0}^{1/2} q_0 \,{{\rm d}}{} \phi}{\theta_0(1/2)+\displaystyle\int_0^{1/2}\left[ \dfrac{\cos{\theta_0}}{F_0^2 q_0^3}-F_0^2q_1-2F_0F_1q_0\right] {{\rm d}}{} \phi + n {\rm \pi}} \quad \text{for} \ n \in \mathbb{Z}^+. \end{equation}

Formula (7.18) provides the crucial eigenvalue condition for the non-existence of solutions. Recall that the ‘parameters’ in this formula, e.g. $\{F_0, F_1, q_0, q_1, \theta _0\}$, are dependent on the chosen energy, $\mathscr {E}$, in (2.9c). Note that $\theta _0(1/2) = 0$ and, in addition, only solutions with positive integer values of $n$ correspond to positive values of the Bond number. Thus, for instance, it is predicted that solutions do not exist at a countably infinite set of discrete values,

(7.19)\begin{equation} B_1(\mathscr{E}) > B_2(\mathscr{E}) > B_3(\mathscr{E}) > \dots > B_n(\mathscr{E}) > \dots > 0. \end{equation}

In the next section, we show that these values of $B$ are associated with points between adjacent ‘fingers’ of solutions in the bifurcation diagram.

Substitution of (7.16) for $C_1$ and $C_{-1}$ into (7.14b) then gives a real-valued solution on the free surface. First, for $\phi <0$, we have

(7.20a)\begin{equation} \bar{\mathfrak{q}}(\phi) \!=\!{-}\frac{2{\rm \pi}}{B^\gamma}\lvert \varLambda_1 \rvert q_0^2 \mathrm{e}^{-{\operatorname{Re}[\chi_1]}/{B}} \left[ \frac{\cos{(G(1/2))}}{\sin (G(1/2))} \cos [\arg{\varLambda_1}\!+\!G(\phi)]\!-\! \sin [\arg{\varLambda_1}\!+\!G(\phi)] \right]\!, \end{equation}

whereas for values on the positive real axis $\phi >0$,

(7.20b)\begin{equation} \bar{\mathfrak{q}}(\phi) \!=\!{-}\frac{2{\rm \pi}}{B^\gamma}\lvert \varLambda_1 \rvert q_0^2 \mathrm{e}^{-{\operatorname{Re}[\chi_1]}/{B}} \left[ \frac{\cos{(G(1/2))}}{\sin (G(1/2))} \cos [\arg{\varLambda_1}\!+\!G(\phi)]\!+\!\sin [\arg{\varLambda_1}\!+\!G(\phi)] \right]\!. \end{equation}

Note that the above forms for $\bar {\mathfrak {q}}$ are valid away from the boundary layer surrounding the Stokes line at $\phi = 0$.

8.1. Finding values for our analytical solution

To obtain precise values for our analytical solution, $\bar {\mathfrak {q}}$, across the domain, we use the form given in (7.14a). This form includes the local change across the boundary layer at $\phi =0$ and requires known values of $q_0$, $\theta _0$, $F_0$, $q_1$ and $F_1$ for a specified value of the energy, $\mathscr {E}$.

In order to calculate values for these nonlinear solutions, we employ Newton iteration on the $O(1)$ and $O(B)$ (4.2) and (4.3) with an even discretisation of the domain, $\phi$. With these values known, the three components of $\bar {\mathfrak {q}}$, the Stokes prefactor, $A_a(\phi )$, the functional prefactor, $Q_a(\phi )$, and the singulant, $\chi _a(\phi )$, may then be calculated individually with a specified value of $B$.

1. (i) For $Q_a(\phi )$ given in (6.7), we take the previously computed values for $\theta _0$, $q_0$, $F_0$, $q_1$ and $F_1$ and employ numerical integration across the domain. As noted in § 6.1, it is convenient to choose a value of $\lvert \varLambda _a \rvert$ in order to facilitate visualisation of the ripples. In figures 9 and 11, we plot $\bar {\mathfrak {q}}$ with $\lvert \varLambda _a \rvert =1$. In figure 10, in order to compare between asymptotic and numerical solutions, we have chosen $\lvert \varLambda _a \rvert =0.006$, which is estimated by numerical fitting. It can be verified that fitting to other fingers changes the constant by only a small amount.

2. (ii) To determine $\chi _a(\phi )$, we split the range of integration as in (6.15). This allows for $\textrm {Re}[\chi _a]$ to be calculated by integrating $q_0$ through the complex-valued domain from the singularity at $f=af^*$ to the wave crest at $f=0$. Next, $\textrm {Im}[\chi _a]$ is found by integrating $q_0$ over the free surface from $f=0$ to $f=\phi$. Values for the integrand, $q_0$, are found with the analytic continuation method from Crew & Trinh (Reference Crew and Trinh2016) described in § 7.1.

3. (iii) To find the Stokes prefactor, $A_a(\phi )$, from (7.12), the upper limit of the integral is determined by using $r_a=\lvert \chi _a \rvert$ and $\bar {\vartheta }={\arg {\chi _a}}/{B}$ from the known values of $\chi _a$. The integral is then calculated with known values of the error function. The constants $C_a$ are then found by calculating $G(1/2)$ from (7.17).

This process yields values for our exponentially small component of the solution, $\bar {\mathfrak {q}}$, for specified values of $B$ and $\mathscr {E}$. The values of $B_n$ where the solvability condition fails from (7.18) are also found with the same method used for $Q_a$ above.

Figure 9.

A comparison between the numerical solution branches of Shelton et al. (Reference Shelton, Milewski and Trinh2021) (shown solid) and the analytical approximations of $B_n$ from (7.18) (shown as black circles). Insets (a)–(d) show the exponentially small ripples, $\bar {\mathfrak {q}}$, from (7.20) for the four locations of $B=0.001876$, $B=0.001264$, $B=0.0009527$ and $B=0.0004978$ (shown as crosses in the main inset). The solutions are all computed at $\mathscr {E} = 0.3804$. A value of $\lvert \varLambda _a \rvert =1$ has been used for the constant prefactor.

Figure 10.

The analytical solution, $q=q_0+Bq_1+\bar {q}$, is shown (line) for the four profiles calculated in figure 9. For comparison, numerical solutions with the same value of $B$ and $\mathscr {E}$ are shown dashed in insets (a) and (b). A value of $\lvert \varLambda _a \rvert = 0.006$ has been used for these comparisons, estimated from numerical comparisons.

8.2. Comparisons

We begin by comparing the values of $B_n$ (where the solvability condition fails) with the ($B,F$) bifurcation space computed numerically by Shelton et al. (Reference Shelton, Milewski and Trinh2021). In taking the same value of the energy, $\mathscr {E}=0.3804$, we visualise these points in the $(B,F)$-plane by approximating $F_n$ by $F_n \approx F_0+B_n F_1$ (an error of $O(B^2)$). This comparison is shown in figure 9. These locations where perturbation solutions are non-existent show excellent agreement with the points between adjacent branches of solutions where numerical solutions could not be calculated.

In addition, four of our analytical solution profiles, $\bar {\mathfrak {q}}$, are shown in insets $(a)$–$(d)$ of this figure. These solutions have been selected to lie in the midpoint of the solution branch, with a Bond number of ($B=(B_n+B_{n+1})/2$). They demonstrate that the ripples obtain their greatest magnitude at the edge of the periodic domain. Note that these ripples are plotted on a zero background state. These same solutions are also shown in figure 10, which includes the first two terms of the asymptotic expansion, $q_0+Bq_1$. These have been provided to compare the magnitude of the ripples in relation to the leading-order Stokes wave.

In our previous numerical work, we demonstrated that as one of the solution branches was transversed, the solution develops an extra wavelength, and this was seen to occur near the top of the solution branch. We observe that the same effect occurs with our analytical solutions. This is demonstrated in figure 11, in which we provide eight solution profiles equally spaced in the Bond number between two adjacent values of $B_n$. From these, we see that as we travel from right to left across the solution branch by decreasing the value of $B$, an additional ripple forms in the centre of the domain.

Figure 11.

Here, for $\mathscr {E}=0.3804$, we plot the exponentially small solution, $\bar {\mathfrak {q}}$, from (7.20) between the two values of $B_{29}=0.0009360$ and $B_{28}=0.0009694$. Note that the base gravity wave is thus not shown. The eight chosen values of $B$ (crosses) are equally spaced between the values of $B_{29}$ and $B_{28}$. This corresponds to the finger $G_{28\to 29}$ found numerically by Shelton et al. (Reference Shelton, Milewski and Trinh2021). A value of $\lvert \varLambda _a\rvert =1$ has been used for the constant prefactor.

8.3. The effects of changing the energy, $\mathscr {E}$

All of the above solutions have been computed for the same fixed value of the energy, $\mathscr {E}=0.3804$. We now relax this restriction by considering values of $\mathscr {E}$ between $0$ and $0.9$. Note that the limiting Stokes wave is not the most energetic (cf. § 6 of Longuet-Higgins & Fox Reference Longuet-Higgins and Fox1978) and for values of $\mathscr {E}$ very close to unity, there are multiple possible solutions beyond the classical Stokes wave. We do not consider solutions too close to the highest wave ($\mathscr {E}>0.9$) in this work.

In figure 12 we show how the locations where the solvability condition fails, $B_n(\mathscr {E})$, change with the energy for values of $n\leqslant 40$. We note that as the energy deceases to zero and we enter the linear regime, these lines tend towards the predictions by Wilton (Reference Wilton1915). These are the discrete values of the Bond number for which two linear solutions of wavenumbers $1$ and $n$ also have the same Froude number. Thus, a single leading-order gravity wave of the type assumed in this work is insufficient for describing Wilton's linear solutions, and is why we recover his values under this limit.

Figure 12.

Values for $B_n$, where the solvability condition fails, are shown for different values of the energy, $\mathscr {E}$. The small-$\mathscr {E}$ predictions by Wilton (Reference Wilton1915) are shown by the black dots at $\mathscr {E}=0$ for $n=20$, $30$ and $40$.

We have also chosen to provide values of $\textrm {Re}[\chi ]$ for different values of $\mathscr {E}$, as this controls the exponential behaviour of the magnitude of our parasitic ripples. This is shown in figure 13, and shows that the constant controlling the exponential behaviour of our solution increases with the energy, $\mathscr {E}$.

Figure 13.

The value of $-\textrm {Re}[\chi ]$ from (6.15) is shown for different values of the energy, $\mathscr {E}$.