Resonant triad interactions of gravity waves in cylindrical basins

Abstract

We present the results of a theoretical investigation into the existence, evolution and excitation of resonant triads of nonlinear free-surface gravity waves confined to a cylinder of finite depth. It is well known that resonant triads are impossible for gravity waves in laterally unbounded domains; we demonstrate, however, that horizontal confinement of the fluid may induce resonant triads for particular fluid depths. For any three correlated wave modes arising in a cylinder of arbitrary cross-section, we prove necessary and sufficient conditions for the existence of a depth at which nonlinear resonance may arise, and show that the resultant critical depth is unique. We enumerate the low-frequency triads for circular cylinders, including a new class of resonances between standing and counter-propagating waves, and also briefly discuss annular and rectangular cylinders. Upon deriving the triad amplitude equations for a finite-depth cylinder of arbitrary cross-section, we deduce that the triad evolution is always periodic, and determine parameters controlling the efficiency of energy exchange. In order to excite a particular triad, we explore the influence of external forcing; in this case, the triad evolution may be periodic, quasi-periodic or chaotic. Finally, our results have potential implications on resonant water waves in man-made and natural basins, such as industrial-scale fluid tanks, harbours and bays.

1. Introduction

Nonlinear resonance is a mechanism by which energy is continuously transferred between a small number of linear wave modes. This phenomenon, first observed in Wilton's analysis of gravity–capillary wave trains (Wilton 1915), has been the subject of frequent investigation over the past century (McGoldrick 1965; Simmons 1969; McGoldrick 1970a,b; Schwartz & Vanden-Broeck 1979; Craik 1986; Hammack & Henderson 1993); indeed, nonlinear resonance has since been observed for wave trains in a growing number of dispersive wave systems, including gravity waves (Phillips 1960; Hasselmann 1961; Benney 1962; Longuet-Higgins 1962), acoustic–gravity waves (Kadri & Stiassnie 2013; Kadri & Akylas 2016), flexural–gravity waves (Wang, Vanden-Broeck & Milewski 2013), two-layer flows (Ball 1964; Joyce 1974; Segur 1980) and atmospheric flows (Raupp et al. 2008; Raupp & Silva Dias 2009). Whilst the aforementioned studies typically consider nonlinear resonance for laterally unbounded domains, the purpose of this study is to demonstrate that energy exchange between free-surface gravity waves may be induced and accentuated by horizontal confinement.

We focus our study on the collective resonance of three linear wave modes, henceforth referred to as a triad (Bretherton 1964). In laterally unbounded domains, the monotonic and concave form of the dispersion curve precludes the existence of resonant triads for gravity wave trains at finite depth (Phillips 1960; Hasselmann 1961), with resonant quartets instead being the smallest possible collective resonant interaction (Benney 1962; Longuet-Higgins 1962; Berger & Milewski 2003). However, confinement of the fluid to a vertical cylinder results in linear wave modes that differ in form to sinusoidal plane waves (except for a rectangular cylinder), so the preclusion of resonant triads no longer applies. Indeed, our study demonstrates that, under certain conditions, resonant triads may arise in cylinders of arbitrary cross-section for specific values of the fluid depth. As resonant triads evolve over a much faster time scale than that of resonant quartets, the exchange of energy in gravity waves is thus more efficient under the influence of lateral confinement (Michel 2019), with potential implications for resonant sloshing in man-made and natural basins (Bryant 1989).

Prior investigations of confined resonant free-surface gravity waves have predominantly focused on the so-called 1:2 resonance, which arises when two of the three linear wave modes comprising a triad coincide. For axisymmetric standing waves in a circular cylinder, Mack (1962) determined a condition for the existence of critical depth-to-radius ratios at which a 1:2 resonance may arise, a result later generalised to cylinders of arbitrary cross-section (Miles 1984b). Miles (1976, 1984a) then characterised the weakly nonlinear evolution of such internal resonances, demonstrating that a 1:2 resonance is impossible in a rectangular cylinder (Miles 1976). Although Miles’ seminal results provide an informative view of the weakly nonlinear dynamics, the influence of fully nonlinear effects was later assessed by Bryant (1989) and Yang et al. (2021). For the case of a circular cylinder of finite depth, Bryant (1989) and Yang et al. (2021) characterised new steadily propagating nonlinear waves arising in the vicinity of a 1:2 resonance, and Yang et al. (2021) also computed nonlinear near-resonant axisymmetric standing waves. Finally, broader mathematical properties of water waves exhibiting O(2) symmetry (of which a circular cylinder is one example) were analysed by Bridges & Dias (1990) and Chossat & Dias (1995).

Given the restrictive set of critical depths at which a 1:2 resonance may arise (Bryant 1989; Yang et al. 2021), it is natural to explore the possibility of nonlinear resonance in cylinders whose depth departs from the depths that trigger a 1:2 resonance. To the best of our knowledge, the first and only such study was the seminal experimental investigation performed by Michel (2019), who focused on resonant triads arising for free-surface gravity waves confined to a finite-depth circular cylinder. Notably, the cylinder depth in Michel's experiment was judiciously chosen so as to isolate a specific triad. Michel utilised bandlimited random horizontal vibration so as to excite two members of the triad, whose nonlinear interaction led to the growth of the third mode. Significantly, the energy of the third mode was, on average, the product of the energies of the remaining two modes, thereby satisfying the quadratic energy exchange typical of resonant triads.

In order to exemplify the mechanism of nonlinear resonance, Michel (2019) also calculated the response of a child mode due to the nonlinear interaction between two parent modes (where all three wave modes comprise the triad). Notably, Michel's calculation is restricted to the early stages of growth and to particular relative phases of the wave modes. In addition, Michel considered a fluid of infinite depth for all but the resonance conditions, for which finite-depth corrections were included. In contrast, we consider general resonances in arbitrary cylinders of finite depth and derive equations for the triad evolution over long time scales. We also believe some nonlinear contributions to the interactions were omitted from Michel's calculation, resulting in quantitative differences (see § 4.3).

The goal of our study is to unify the existence and evolution of 1:2 and triadic resonances into a single mathematical framework, effectively characterising all triad interactions of this type. Based on existing theory, it is unclear how the existence of resonant triads depends on the form of the cylinder cross-section, and which combinations of wave modes are permissible for judicious choice of the fluid depth. Furthermore, the range of depths that may excite a particular triad is uncertain, with 1:2 resonances only excited in a very narrow window about each critical depth (Mack 1962; Miles 1984b). Once a particular triad is excited, one anticipates that the triad evolution will be governed by the canonical triad equations (Bretherton 1964; Craik 1986); however, quantifying the triad evolution and relative energy exchange requires computation of the triad coupling coefficients. Finally, it is unclear how best to excite triads in arbitrary cylinders, both with and without external forcing.

We here present a relatively comprehensive characterisation of the existence, evolution and excitation of resonant triads for gravity waves confined to a cylinder of arbitrary cross-section and finite depth. In order to reduce the problem to its key components, we first truncate the Euler equations, recasting the fluid evolution in terms of a finite-depth Benney–Luke equation (§ 2), incorporating only the nonlinear interactions necessary for resonant triads. In § 3, we prove necessary and sufficient conditions for there to exist a finite depth at which three linear wave modes may form a resonant triad. In particular, we prove that resonant triads are impossible for rectangular cylinders, yet there is an abundance of resonant triads for circular cylinders. We then use multiple-scale analysis to determine the long-time evolution of a triad in a cylinder of arbitrary cross-section (§ 4), from which we characterise the relative coupling of different triads. Finally, we explore the excitation of resonant triads (§ 5), and discuss the potential extension of our theoretical developments to the cases of applied forcing and two-layer flows (§ 6).

2. Formulation

We consider the irrotational flow of an inviscid, incompressible liquid that is bounded above by a free surface, confined laterally by the vertical walls of a cylinder whose horizontal cross-section, $\mathcal {D}$, is enclosed by the curve $\partial \mathcal {D}$ and bounded below by a rigid horizontal plane lying a distance $H$ below the undisturbed free surface; see figure 1. We consider the fluid evolution in dimensionless variables, taking the cylinder's typical horizontal extent, $a$, as the unit of length, and $\sqrt {ag^{-1}}$ as the unit of time, where $g$ is the acceleration due to gravity. It follows that the dimensionless free-surface elevation, $\eta (\boldsymbol {x},t)$, and velocity potential, $\phi (\boldsymbol {x},z,t)$, evolve according to the equations

(2.1a)$$\begin{gather} \Delta \phi + \phi_{zz} = 0 \quad \mathrm{for}\ \boldsymbol{x} \in \mathcal{D},\quad -h < z < \epsilon \eta, \end{gather}$$

(2.1b)$$\begin{gather}\phi_t + \eta + \frac{\epsilon}{2}(|\boldsymbol{\nabla} \phi|^2 + \phi_z^2) = 0 \quad \mathrm{for}\ \boldsymbol{x} \in \mathcal{D}, \quad z = \epsilon \eta, \end{gather}$$

(2.1c)$$\begin{gather}\eta_t + \epsilon \boldsymbol{\nabla}\phi \boldsymbol{\cdot} \boldsymbol{\nabla}\eta = \phi_z \quad \mathrm{for}\ \boldsymbol{x} \in \mathcal{D}, \quad z = \epsilon \eta, \end{gather}$$

(2.1d)$$\begin{gather}\boldsymbol{n}\boldsymbol{\cdot}\boldsymbol{\nabla} \phi = 0 \quad \mathrm{for}\ \boldsymbol{x} \in \partial\mathcal{D}, \quad -h < z < \epsilon \eta, \end{gather}$$

(2.1e)$$\begin{gather}\phi_z = 0 \quad \mathrm{for}\ \boldsymbol{x} \in \mathcal{D},\quad z ={-}h, \end{gather}$$

corresponding to the continuity equation, dynamic and kinematic boundary conditions and no flux through the vertical walls and horizontal base, respectively. In (2.1), the dimensionless parameter $\epsilon$ is proportional to the typical wave slope, $h = H/a$ is the ratio of the fluid depth to the typical horizontal extent, $\boldsymbol {n}$ is a unit vector normal to the boundary $\partial \mathcal {D}$ and the operators $\boldsymbol {\nabla }$ and $\varDelta$ denote the horizontal gradient and Laplacian, respectively. Moreover, conservation of mass implies that the free surface satisfies $\iint _{\mathcal {D}} \eta \,\mathrm {d}A = 0$ for all time. Finally, in dimensional variables, $a\boldsymbol {x}$ is the two-dimensional horizontal coordinate, $az$ is the upward-pointing vertical coordinate, $\sqrt {ag^{-1}} t$ denotes time, $\epsilon a \eta$ is the free-surface displacement and $\epsilon a \sqrt {ag}\phi$ is the velocity potential.

Figure 1. Schematic diagram of the cylindrical tank (with cross-section $\mathcal {D}$ and boundary $\partial \mathcal {D}$) partially filled with liquid. The undisturbed free surface (dashed lines) lies on $z = 0$, a distance $H$ above the rigid bottom plane (grey). The disturbed free surface is sketched in dash-dotted lines.

We aim to develop a broad framework for understanding resonant triads in a cylinder of finite depth; however, care must be taken when modelling fluid–boundary interactions and determining the class of permissible cylinder cross-sections. From a modelling perspective, we employ an assumption generally implicit to the water-wave problem in bounded domains; specifically, we neglect the meniscus and dissipation arising near the vertical walls (Miles 1967), thus determining that the free surface intersects the boundary normally, i.e. $\boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {\nabla } \eta = 0$ for $\boldsymbol {x} \in \partial \mathcal {D}$ (Miles & Henderson 1990). In order to maximise the generality of our investigation, we allow the cylinder cross-section, $\mathcal {D}$, to be fairly arbitrary; however, the mathematical developments presented herein require $\mathcal {D}$ to be bounded with a piecewise-smooth boundary, thereby allowing us to utilise the spectral theorem for compact self-adjoint operators (Kreyszig 1989) and the divergence theorem. As most cylinders of practical interest consist of a piecewise-smooth boundary, this mathematical restriction fails to limit the breadth of our study.

2.1. Derivation of the Benney–Luke equation

As our study is focused on the weakly nonlinear evolution of small-amplitude waves, we proceed to simplify (2.1) in the case $0 < \epsilon \ll 1$ and $h = O(1)$. We begin by expanding the dynamic and kinematic boundary conditions (2.1b)–(2.1c) about $z = 0$ in powers of $\epsilon$, which, upon eliminating $\eta$, gives rise to the equation (Benney 1962; Milewski & Keller 1996)

(2.2)\begin{align} \phi_{tt} + \phi_z = \epsilon(\partial_t (\phi_{tz}\phi_t) + \phi_{zz}\phi_t - \partial_t(|\boldsymbol{\nabla} \phi|^2) - \phi_z\phi_{zt}) + O(\epsilon^2) \quad\mathrm{for}\ \boldsymbol{x} \in \mathcal{D}, \quad z = 0. \end{align}

To reduce the fluid evolution to the dynamics arising on the linearised free surface, $z = 0$, we define the Dirichlet-to-Neumann operator, $\mathscr {L}$, so that $\mathscr {L} \phi |_{z = 0} = \phi _z|_{z = 0}$. Here, $\phi$ satisfies Laplace's equation (2.1a) over the linearised domain $-h < z < 0$, with $\partial _z \phi = 0$ on $z = -h$ (see (2.1e)) and $\boldsymbol {n}\boldsymbol {\cdot } \boldsymbol {\nabla } \phi = 0$ for $\boldsymbol {x} \in \partial \mathcal {D}$ (see (2.1d)). Notably, the Dirichlet-to-Neumann operator may be defined in terms of its spectral representation, as detailed in § 2.2. By denoting $u(\boldsymbol {x},t) = \phi (\boldsymbol {x},0,t)$, we finally obtain the finite-depth Benney–Luke equation (Benney 1962; Benney & Luke 1964; Milewski & Keller 1996)

(2.3)\begin{equation} u_{tt} + \mathscr{L} u + \epsilon\left( u_t(\mathscr{L}^2 + \varDelta)u + \frac{\partial {}}{\partial {t}} [(\mathscr{L} u)^2 +|\boldsymbol{\nabla} u|^2]\right) = O(\epsilon^2) \quad\mathrm{for}\ \boldsymbol{x} \in \mathcal{D}, \end{equation}

where we have simplified the nonlinear terms in (2.2) using $\phi _{zz} = -\Delta \phi$ and $u_{tt} = -\mathscr {L} u + O(\epsilon )$.

The remainder of our investigation will be focused on the evolution of resonant triads governed by the Benney–Luke equation (2.3). As resonant triads arising in confined geometries are governed primarily by quadratic nonlinearities, it is sufficient to neglect terms of size $O(\epsilon ^2)$ in (2.3); however, higher-order corrections to the Benney–Luke equation may be derived by following a similar expansion procedure (Benney 1962; Milewski & Keller 1996; Berger & Milewski 2003). Although our investigation is mainly focused on the evolution of the velocity potential, $u$, one may recover the leading-order free-surface elevation from the dynamic boundary condition (2.1b), namely $\eta = -u_t + O(\epsilon )$.

2.2. Spectral representation of the Dirichlet-to-Neumann operator

The Dirichlet-to-Neumann operator, $\mathscr {L}$, may be understood in terms of the discrete set of orthogonal eigenfunctions of the horizontal Laplacian operator (Kreyszig 1989). Specifically, we consider the set of real-valued eigenfunctions, $\varPhi _n(\boldsymbol {x})$, satisfying

(2.4)\begin{equation} -\Delta \varPhi_n = k_n^2 \varPhi_n \quad \mathrm{for}\ n = 0, 1, \ldots, \end{equation}

where the corresponding eigenvalues, $k_n^2$, are ordered so that $0 = k_0 < k_1\leq k_2 \leq \ldots$. Moreover, each eigenfunction satisfies the boundary condition $\boldsymbol {n}\boldsymbol {\cdot }\boldsymbol {\nabla } \varPhi _n = 0$ on $\partial \mathcal {D}$, as motivated by the no-flux condition (2.1d). Finally, the orthogonal eigenfunctions are normalised so that $\langle \varPhi _m, \varPhi _n \rangle = \delta _{mn}$, where

(2.5)\begin{equation} \langle \,f, g\rangle = \frac{1}{S}\iint_\mathcal{D} f g \,\mathrm{d}A \end{equation}

defines an inner product for real functions $f$ and $g$, $S$ is the area of $\mathcal {D}$ and $\delta _{mn}$ is the Kronecker delta. Notably, $\varPhi _0(\boldsymbol {x}) = 1$ is the constant eigenfunction, with corresponding eigenvalue $k_0 = 0$.

To determine the Dirichlet-to-Neumann operator for sufficiently smooth $\phi$, we first substitute the series expansion $\phi (\boldsymbol {x},z) = \sum _{n = 0}^\infty \phi _n(z)\varPhi _n(\boldsymbol {x})$ into Laplace's equation (2.1a), where we have temporally omitted the time dependence. We then solve the resulting equation for $\phi _n(z)$ over the linearised domain $-h < z < 0$, in conjunction with the no-flux condition on $z = -h$ (see (2.1e)). It follows that $\partial _z\phi _n(0) = \hat {\mathscr {L}}_n \phi _n(0)$, where

(2.6)\begin{equation} \hat{\mathscr{L}}_n = k_n \tanh(k_n h) \end{equation}

is the spectral multiplier of the Dirichlet-to-Neumann operator, $\mathscr {L}$. By expressing the time-dependent free-surface velocity potential, $u = \phi |_{z = 0}$, in terms of the basis expansion $u(\boldsymbol {x}, t ) = \sum _{n = 0}^\infty u_n(t) \varPhi _n(\boldsymbol {x})$, it follows that the Dirichlet-to-Neumann map has the spectral representation $\mathscr {L} u = \sum _{n = 0}^\infty \hat {\mathscr {L}}_n u_n \varPhi _n$.

3. The existence of resonant triads

Resonant triads arise due to the exchange of energy between linear wave modes, an effect induced by nonlinear wave interactions. In order to define resonant triads mathematically, it is necessary to first determine the angular frequency associated with each linear wave mode. In the limit $\epsilon \rightarrow 0$, the Benney–Luke equation (2.3) reduces to the linear equation $u_{tt} + \mathscr {L} u = 0$. By seeking a solution to the linearised Benney–Luke equation of the form $u(\boldsymbol {x},t) = \varPhi _n(\boldsymbol {x}) \exp ({-\mathrm {i} \omega _n t})$, we conclude that the angular frequency, $\omega _n$, satisfies $\omega _n^2 = \hat {\mathscr {L}}_n$, or the more familiar (Lamb 1932)

(3.1)\begin{equation} \omega_n^2 = k_n \tanh(k_n h). \end{equation}

As we will see, a crucial aspect of the following analysis is that the angular frequency depends on the fluid depth, i.e. $\omega _n(h)$. Finally, we note that the angular frequency is larger for more oscillatory eigenfunctions (i.e. for larger values of $k_n$); by analogy to the evolution of plane gravity waves, we refer to $k_n$ as a ‘wavenumber’ henceforth.

We proceed by considering three linear wave modes, enumerated $n_1$, $n_2$ and $n_3$, where we denote

(3.2a–c)\begin{equation} \varOmega_j = \omega_{n_j}, \quad K_j = k_{n_j}, \quad \mathrm{and}\quad \varPsi_j(\boldsymbol{x}) = \varPhi_{n_j}(\boldsymbol{x}) \quad \mathrm{for}\ j = 1, 2, 3. \end{equation}

Notably, we exclude the wavenumber $k_0 = 0$ from consideration as the corresponding eigenmode, $\varPhi _0$, simply reflects the invariance of the Benney–Luke equation (2.3) under the mapping $u \mapsto u + \mathrm {constant}$; henceforth, we consider only wavenumbers $K_j > 0$. The three linear wave modes form a resonant triad if there is a critical fluid depth, $h_c$, satisfying

(3.3)\begin{equation} \varOmega_1(h_c) \pm \varOmega_2(h_c) \pm \varOmega_3(h_c) = 0, \end{equation}

where all four sign combinations are permissible (we consider $\varOmega _j > 0$ without loss of generality). To simplify notation in the following arguments, we restrict our attention to the particular case

(3.4)\begin{equation} \varOmega_1(h_c) + \varOmega_2(h_c) = \varOmega_3(h_c), \end{equation}

where the other three sign combinations in (3.3) may be recovered by suitable re-indexing of the $\varOmega _j$ terms. However, as we will see in § 4, an additional constraint necessary for triads to exist is the eigenmode correlation condition,

(3.5)\begin{equation} \iint_{\mathcal{D}} \varPsi_1\varPsi_2\varPsi_3 \,\mathrm{d}A \neq 0, \end{equation}

which implies that the product of any two eigenmodes is non-orthogonal to the remaining eigenmode.

3.1. The existence of a critical depth

We proceed to determine necessary and sufficient conditions on the wavenumbers, $K_j$, for there to exist a depth, $h_c$, at which a resonant triad forms, where such a critical depth is unique. We summarise our results in terms of the following theorem:

Theorem 1 There exists a positive and finite value of $h$ such that $\varOmega _1 + \varOmega _2 = \varOmega _3$ if and only if

(3.6)\begin{equation} K_1 + K_2 < K_3 < (\sqrt{K_1} + \sqrt{K_2})^2. \end{equation}

When this pair of inequalities is satisfied, the corresponding value of $h$ is unique.

We briefly sketch the proof of Theorem 1, with full details presented in Appendix A. We first demonstrate that no solutions to $\varOmega _1 + \varOmega _2 = \varOmega _3$ are possible when the bounds in (3.6) are violated, i.e. when $K_1 + K_2 \geq K_3$ or when $\sqrt {K_1} + \sqrt {K_2} \leq \sqrt {K_3}$. We then consider the case where the inequalities (3.6) are satisfied and determine the existence of positive roots to the function $F(h) = (\varOmega _1(h) +\varOmega _2(h))/\varOmega _3(h) - 1$. In this case, we demonstrate that $\lim _{h \rightarrow 0} F(h) < 0$ and $\lim _{h\rightarrow \infty } F(h) > 0$, from which we conclude that $F(h)$ has at least one root (by continuity of $F$). Finally, we deduce that this root is unique by proving that $F(h)$ is a strictly monotonically increasing function of $h$ when the inequalities (3.6) are satisfied.

Two important conclusions may be deduced from Theorem 1. First, it follows from (3.6) that the wavenumber, $K_3$, corresponding to the largest angular frequency, $\varOmega _3$, is larger than both the other two wavenumbers ($K_1$ and $K_2$), but it cannot be arbitrarily large (as supplied by the upper bound). For a given pair of eigenmodes (say $\varPsi _1$ and $\varPsi _2$), we conclude that there are likely to be only finitely many eigenmodes that can resonate with this pair (indeed, that number might fairly small, or even zero). Second, when modes 1 and 2 coincide (a 1:2 resonance), one deduces that $\varOmega _1 = \varOmega _2$ and $K_1 = K_2$; as such, the existence bounds (3.6) simplify to $2K_1 < K_3 < 4K_1$, or $2 < K_3/K_1 < 4$ (Mack 1962; Miles 1984b).

3.2. Determining the critical depth

Although Theorem 1 determines necessary and sufficient conditions on the wavenumbers, $K_j$, for there to be a critical depth, $h_c$, at which a resonant triad exists, the critical depth remains to be determined. In general, the critical depth must be computed numerically (being the unique root of the nonlinear function $F(h)$); however, we demonstrate that useful quantitative and qualitative information may be obtained via asymptotic analysis. For the remainder of this section, we consider the rescaled wavenumbers, $\xi _1 = K_1/K_3$ and $\xi _2 = K_2/K_3$, and the rescaled depth, $\zeta = K_3 h$; it remains to determine the root, $\zeta _c$, of

(3.7)\begin{equation} F(\zeta) = \sqrt{\frac{ \xi_1 \tanh(\xi_1 \zeta )} {\tanh(\zeta)} } + \sqrt{\frac{ \xi_2 \tanh(\xi_2 \zeta )} {\tanh(\zeta)} } - 1, \end{equation}

when $\xi _1, \xi _2 > 0$ satisfy

(3.8)\begin{equation} \xi_1 + \xi_2 < 1 < \sqrt{\xi_1} + \sqrt{\xi_2}. \end{equation}

In figure 2(a), we present contours of the critical rescaled depth, $\zeta _c$, in the $(\xi _1, \xi _2)$-plane, restricted to the region demarcated by (3.8). Consistent with the limits $\lim _{\zeta \rightarrow 0}F(\zeta ) = \xi _1 + \xi _2 - 1$ and $\lim _{\zeta \rightarrow \infty }F(\zeta ) = \sqrt {\xi _1} + \sqrt {\xi _2} - 1$, we observe that the root, $\zeta _c$, tends to zero at the line $\xi _1 + \xi _2 = 1$, and approaches infinity at the curve $\sqrt {\xi _1} + \sqrt {\xi _2} = 1$. Furthermore, the uniqueness of the root of $F$ for given $(\xi _1, \xi _2)$ is reflected in the observation that the contours of $\zeta _c$ do not cross. Finally, we note that the contours are symmetric about the line $\xi _1 = \xi _2$, which is a direct consequence of the invariance of $F(\zeta )$ under the mapping $\xi _1 \leftrightarrow \xi _2$ (see (3.7)).

Figure 2. Contours of the rescaled critical depth, $\zeta _c = h_cK_3$, as a function of the rescaled wavenumbers, $\xi _1 = K_1/K_3$ and $\xi _2 = K_2/K_3$. (a) The contours computed numerically from (3.7). The black lines indicate the limiting cases of $\zeta _c \rightarrow 0$ (at $\xi _1 + \xi _2 = 1$) and $\zeta _c \rightarrow \infty$ (at $\sqrt {\xi _1} + \sqrt {\xi _2} = 1$). (b) The contours are overlaid by the leading-order approximation ((3.15); circles) and the higher-order correction ((3.16); diamonds) for $\zeta _c$ equal to 0.5, 1, 1.5 and 2.

Although we are primarily interested in the physically relevant case for which the cylinder's depth-to-width ratio, $h$, is of size $O(1)$, an informative analytic result may be obtained by considering $F(\zeta )$ in the limit $\zeta \ll 1$ (or $K_3 h \ll 1$). By utilising the Taylor expansion

(3.9)\begin{equation} \sqrt{\tanh(x)} \sim \sqrt{x}\left(1 - \frac{x^2}{6} + \frac{19}{360}x^4 + O(x^6)\right), \end{equation}

we obtain

(3.10)\begin{align} & \sqrt{\xi_1\tanh(\xi_1\zeta)} + \sqrt{\xi_2\tanh(\xi_2\zeta)} -\sqrt{\tanh(\zeta)} \nonumber\\ &\quad \sim \sqrt{\zeta}\left[\left(\xi_1 + \xi_2 - 1\right) - \frac{\zeta^2}{6}(\xi_1^3 + \xi_2^3 - 1) + O(\zeta^4)\right] , \end{align}

for $0 < \zeta \ll 1$. Whilst deriving (3.10), we have utilised the bound $\xi _1, \xi _2 < 1$ (see (3.8)), which additionally ensures that $0 < \xi _j\zeta \ll 1$ for $j = 1,2$. We note that the left-hand side of (3.10) is equal to $F(\zeta )\tanh (\zeta )$, so $\zeta _c$ satisfies

(3.11)\begin{equation} \xi_1 + \xi_2 - 1 - \frac{\zeta_c^2}{6}(\xi_1^3 + \xi_2^3 - 1) = O(\zeta_c^4), \end{equation}

provided that $0 < \zeta _c \ll 1$. By neglecting terms of size $O(\zeta _c^4)$ in (3.11), one may then easily solve for $\zeta _c$ in terms of $\xi _1$ and $\xi _2$.

Alternatively, a more succinct expression for $\zeta _c$ may be found by first noting that

(3.12)\begin{equation} \xi_1^3 + \xi_2^3 = (\xi_1 + \xi_2)^3 - 3\xi_1\xi_2(\xi_1 + \xi_2) = 1 - 3\xi_1\xi_2 + O(\zeta_c^2), \end{equation}

where we have utilised the leading-order approximation $\xi _1 + \xi _2 = 1 + O(\zeta _c^2)$ (see (3.11)) to determine the second equality. Upon substituting (3.12) into (3.11), we find that $\xi _1$, $\xi _2$ and $\zeta _c$ are now related by the notably simpler expression

(3.13)\begin{equation} \xi_1 + \xi_2 - 1 + \frac{\zeta_c^2}{2}\xi_1\xi_2 = O(\zeta_c^4). \end{equation}

By neglecting terms of $O(\zeta _c^4)$, the leading-order approximation for the rescaled critical depth, $\zeta _c$, is given by

(3.14)\begin{equation} \zeta_c \sim \sqrt{\frac{2(1 - \xi_1 - \xi_2)}{\xi_1\xi_2}}, \end{equation}

an expression valid when $0 < \zeta _c \ll 1$ and $\xi _1 + \xi _2 < 1$ (see (3.8)). Alternatively, one may deduce from (3.13) that the contours of $\zeta _c$ satisfy the approximate form

(3.15)\begin{equation} \xi_2 \sim \frac{1 - \xi_1}{1 + \dfrac{1}{2}\zeta_c^2\xi_1}, \end{equation}

where the term in the denominator is responsible for the increased ‘bending’ of the contours as $\zeta _c$ becomes progressively larger (see figure 2). We note that the additional simplification afforded by (3.12) allows for a far more tractable representation of the contours relative to solving (3.11) directly for $\xi _2$ given $\xi _1$ and $\zeta _c$.

Despite being derived under the assumption $0 < \zeta _c \ll 1$, we see in figure 2(b) that the contours given by (3.15) agree favourably with the numerical solution even up to $\zeta _c \approx 1$. However, it is readily verified from (3.14) that the asymptotic approximation of each contour crosses the boundary curve $\sqrt {\xi _1} + \sqrt {\xi _2} = 1$ at $\zeta _c = 4$ (for which $\xi _1 = \xi _2 = \frac {1}{4}$), thereby demonstrating that the reduced asymptotic form has limited applicability (even in a qualitative sense) for slightly larger values of $\zeta _c$. One may further improve the quantitative (and, to an extent, qualitative) agreement between the asymptotic analysis and numerical computation by including terms of size $O(\zeta ^4)$ in (3.10); indeed, an analogous calculation gives rise to the following higher-order correction to (3.13):

(3.16)\begin{equation} \xi_1 + \xi_2 - 1 + \frac{\zeta_c^2}{2}\xi_1\xi_2 + \frac{\zeta_c^4}{72} \xi_1\xi_2\left(\xi_1\xi_2 - 1\right) = O(\zeta_c^6). \end{equation}

Although one may then solve for $\zeta _c$ given $\xi _1$ and $\xi _2$ (or, alternatively, determine the contours of $\zeta _c$) by truncating terms of $O(\zeta _c^6)$ in (3.16), the resulting algebraic expressions yield little qualitative information. However, one may, in principle, use this reduced form as a reasonable initial guess for a numerical root-finding algorithm for determining the root of $F(\zeta )$, provided that $\zeta _c$ is not too large.

3.3. Example cavities

Our investigation into the emergence of resonant triads has been focused, thus far, on finite-depth cylinders with arbitrary horizontal cross-section. However, it is convenient to understand how the results of Theorem 1 influence the formation (or not) of resonant triads for some specific cross-sections, namely rectangular, circular and annular cylinders.

3.3.1. Rectangular cylinder

It is well known that resonant triads are impossible for plane gravity waves evolving across an unbounded horizontal domain of finite depth (Phillips 1960; Hasselmann 1961) (weak interactions are possible, however, in the shallow-water limit, $K_jh\rightarrow 0$, for which $\tanh (K_jh)$ in the dispersion relation (3.1) is replaced by its leading-order approximation, $K_j h$ (Phillips 1960; Bryant 1973; Miles 1976)). We now utilise Theorem 1 to demonstrate a similar result: resonant triads are impossible for gravity waves evolving within a rectangular cylinder of finite depth. Our result generalises the special case of a 1:2 resonance, for which the impossibility of internal resonance in a rectangular cylinder was demonstrated by Miles (1976).

To proceed, we consider a rectangular cylinder with side lengths $L_x$ and $L_y$. By orientating the Cartesian coordinate system, $\boldsymbol {x} = (x,y)$, so that the cylinder cross-section is defined by the region $0 < x < L_x$ and $0 < y < L_y$, the eigenmodes are of the form

(3.17)\begin{equation} \varPhi_{mn}(x,y) = \frac{1}{\mathscr{N}_{mn}}\cos(p_m x)\cos(q_n y), \end{equation}

where $\mathscr {N}_{mn} > 0$ is a normalisation constant. Notably, the wavenumbers $p_m = m{\rm \pi} /L_x$ and $q_n = n{\rm \pi} /L_y$ are chosen so that the no-flux condition is satisfied (see (2.1d)). For a triad determined by the non-negative integers $m_j$ and $n_j$ (for $j = 1, 2, 3$), the corresponding wavenumbers, $P_j = p_{m_j}$ and $Q_j = q_{n_j}$, must satisfy $P_1 + P_2 = P_3$ and $Q_1 + Q_2 = Q_3$ (under suitable reordering of the subscripts) in order for the eigenmode correlation condition (3.5) to be satisfied. By defining the wave vector $\boldsymbol {k}_j = (P_j, Q_j)$, the conditions on $P_j$ and $Q_j$ simplify to the single requirement $\boldsymbol {k}_1 + \boldsymbol {k}_2 = \boldsymbol {k}_3$, where the triangle inequality supplies that $|\boldsymbol {k}_3| \leq |\boldsymbol {k}_1| + |\boldsymbol {k}_2|$. As the eigenvalues, $K_j^2$, of the negative Laplacian operator are related to the wave vectors via $K_j = |\boldsymbol {k}_j|$, we deduce that $K_3 \leq K_1 + K_2$. Owing to the violation of the left-hand bound in (3.6), we conclude that resonant triads cannot exist in a rectangular cylinder of finite depth.

3.3.2. Circular cylinder

We consider a circular cylinder of unit radius in dimensionless variables (i.e. the dimensional radius is equal to $a$; see § 2). For polar coordinates $\boldsymbol {x} = (r,\theta )$, it is well known that the corresponding (complex-valued) eigenmodes may be expressed in the form

(3.18)\begin{equation} \varPhi_{mn}(r,\theta) = \frac{1}{\mathscr{N}_{mn}}\mathrm{J}_m(k_{mn}r)\mathrm{e}^{\mathrm{i} m \theta},\quad\mathrm{where}\ \mathscr{N}_{mn} = |\mathrm{J}_m(k_{mn})|\sqrt{1 - \frac{m^2}{k_{mn}^2}} \end{equation}

is the normalisation factor and $m$ is the azimuthal wavenumber (an integer). Furthermore, the no-flux condition (2.1d) determines that the radial wavenumbers, denoted $k_{mn}$, satisfy $\mathrm {J}_m'(k_{mn}) = 0$, where $0 < k_{m1} < k_{m2} < \cdots$ (we exclude $k_{00} = 0$ from consideration; see § 3). Notably, the eigenvalues of the negative Laplacian operator are precisely the squared wavenumbers, $k_{mn}^2$; consequently, the antinodes of each Bessel function play a pivotal role in determining the existence of resonant triads.

Akin to the rectangular cylinder, we find that the eigenmode correlation condition imparts an important restriction on the combination of eigenmodes that may resonate. For given $m_j$ and $n_j$ (for $j = 1,2,3$), we denote $K_j = k_{m_jn_j}$, $\varPsi _j = \varPhi _{m_j n_j}$ and $N_j = \mathscr {N}_{m_j n_j}$. Although the correlation condition given in (3.5) is defined for real eigenmodes, a similar condition holds for complex-valued eigenmodes, namely $\iint _{\mathcal {D}} \varPsi _1 \varPsi _2 \varPsi _3^*\,\mathrm {d}A \neq 0$. By considering the quantity

(3.19)\begin{align} \iint_{\mathcal{D}} \varPsi_1 \varPsi_2 \varPsi_3^* \,\mathrm{d}A &= \frac{1}{N_1 N_2 N_3}\left(\int_0^1 r\mathrm{J}_{m_1}(K_1 r) \mathrm{J}_{m_2}(K_2 r) \mathrm{J}_{m_3}(K_3 r)\,\mathrm{d}r\right) \nonumber\\ &\quad \times \left(\int_0^{2{\rm \pi}} \exp({\mathrm{i}(m_1 + m_2 - m_3)\theta})\,\mathrm{d}\theta\right), \end{align}

we deduce from the azimuthal integral that a necessary condition for the correlation integral to be non-zero is $m_1 + m_2 = m_3$ (Michel 2019). This condition thus restricts the permissible combinations of azimuthal wavenumbers in a manner similar to the restriction on the permissible planar wavenumbers for the case of a rectangular cylinder. Unlike rectangular cylinders, however, we demonstrate that resonant triads are possible in a circular cylinder.

Despite the apparent restriction of the Bessel antinodes, $K_j$, and summation condition on the azimuthal wavenumbers, $m_j$, Theorem 1 determines that a vast array of resonant triads may be excited for judicious choices of the fluid depth. In table 1, we list a small number of resonant triads and each corresponding critical depth, $h_c$, subject to the restrictions $|m_j| \leq 3$ and $n_j \leq 3$; for larger values of $|m_j|$ and $n_j$, the corresponding wave field becomes increasingly oscillatory, to the extent that the effects of surface tension and dissipation might become appreciable. Moreover, even marginally relaxing the upper bounds on $|m_j|$ and $n_j$ vastly increases the number of resonant triads; indeed, the restriction $|m_j| \leq 4$ and $n_j \leq 4$ introduces 70 additional resonant triads relative to table 1. As the upper bounds for $|m_j|$ and $n_j$ are further increased, the typical difference between the various critical depths decreases and an increasingly large number of triads form at small values of the critical depth. Triads forming in shallow fluids (e.g. triads 21 to 25 in table 1) have physical relevance only at larger length scales (e.g. lakes) as dissipation could become a dominant factor at smaller scales.

Table 1. Combinations of the azimuthal wavenumbers, $m_j$, and radial mode indices, $n_j$, that form a resonant triad ($m_1 + m_2 = m_3$ and $\varOmega _1 + \varOmega _2 = \varOmega _3$) at critical depth, $h_c$, in a circular cylinder of unit radius. For each triad, the corresponding wavenumbers, $K_j = k_{m_j n_j}$, satisfy (3.6), and the correlation condition, $\iint _{\mathcal {D}} \varPsi _1\varPsi _2\varPsi _3^*\,\mathrm {d}A \neq 0$, is met. The list is restricted to resonant triads arising for $|m_j|, n_j \leq 3$, and we consider $m_1 \leq m_2$ and $m_3 \geq 0$ without loss of generality. We have omitted resonances that give rise to the same critical depth, but with the roles of modes 1 and 2 swapped. The triad numbers (left column) and shaded rows are referenced in the text.

Although the list of triads in table 1 is restricted to the lowest radial and azimuthal modes, we observe some general trends. In particular, we observe that the correlation integral, $\iint _{\mathcal {D}} \varPsi _1\varPsi _2\varPsi _3^*\,\mathrm {d}A$, generally decreases in magnitude as the fluid depth increases. Although the correlation integral remains non-zero (as is necessary to satisfy the correlation condition), its small value in some cases (e.g. triad 10) potentially corresponds to an elongation of the triad evolution time scale (see § 4.1). Moreover, we observe that the average wavenumber involved in the triad, $\bar {K} = \frac {1}{3}(K_1 + K_2 + K_3)$, is appreciably larger when the critical depth is very small. This correlation is consistent with the form of the corresponding angular frequency, $\varOmega _j = \sqrt {K_j \tanh (K_j h_c)}$, for which a small critical depth, $h_c$, is necessary for finite-depth effects to be appreciable when the typical wavenumber is large. To enumerate the myriad resonant triads arising in a circular cylinder when the upper bounds on $|m_j|$ and $n_j$ are relaxed, we provide MATLAB code in the supplementary material available at https://doi.org/10.1017/jfm.2023.441.

At this juncture, it is informative to assess how the triads listed in table 1 relate to the resonances explored in prior investigations. First, triad 8 in table 1 (dark grey row) was explored by Michel (2019) for a circular cylinder of radius 9.45 cm and an approximate fluid depth of 3 cm; it follows that the depth-to-radius ratio in Michel's experiment was approximately 0.317, close to the value of 0.30197 reported in table 1. Furthermore, table 1 (grey rows) incorporates two well-known examples of a 1:2 resonance, for which modes 1 and 2 coincide: (i) the critical depth $h_c = 0.83138$ (triad 2) corresponds to the second-harmonic resonance with the fundamental mode (Miles 1976, 1984a; Bryant 1989; Yang et al. 2021); (ii) the critical depth $h_c = 0.19814$ (triad 14) corresponds to a standing wave composed of two resonant axisymmetric modes (Mack 1962; Yang et al. 2021). Finally, triads 1, 6, 16, 19 and 20 (table 1, light grey rows) form an interesting class of resonant triad, for which an axisymmetric mode ($m_3 = 0$) interacts with two identical counter-propagating non-axisymmetric modes ($m_1= -m_2 \neq 0$ and $n_1 = n_2$). In fact, our investigation in § 5.1 demonstrates that the axisymmetric mode is the so-called pump mode, and may thus excite the non-axisymmetric modes, even when the initial energy in each non-axisymmetric mode is negligible. We draw an analogy between this novel class of resonant triad and the excitation of beach edge waves (Guza & Davis 1974) in § 6.

We conclude our exploration of resonant triads arising in a circular cylinder by remarking that the fluid depth may, in some cases, be judiciously chosen so as to excite multiple triads. In general, the condition on the angular frequencies, $\varOmega _1 + \varOmega _2 = \varOmega _3$ (see (3.4)), cannot be satisfied for two distinct triads at the same fluid depth; however, nonlinear resonance may persist for both triads provided that each condition on the angular frequencies is approximately satisfied (Bretherton 1964; McGoldrick 1972; Craik 1986), at the cost of weak detuning (see § 4.3.1 for further details). Specifically, if triads 1 and 2 have critical depths $h_{c,1}$ and $h_{c,2}$, respectively, then there is potential excitement of both triads when the fluid depth, $h$, satisfies $|h - h_{c,j}| = O(\epsilon )$ for $j = 1,2$ (where $0 < \epsilon \ll 1$ is the typical wave slope; see § 2), giving rise to the approximation $\varOmega _1 + \varOmega _2 - \varOmega _3 = O(\epsilon )$ for each triad. For example, if $0 < h_{c,2} - h_{c,1} \ll 1$, then it may be sufficient to excite both triads at an intermediate depth, $h_{c,1} \leq h \leq h_{c,2}$. We note, however, that the excitation of multiple triads at a single fluid depth is not possible when the depth discrepancy, $|h - h_{c,j}|$, becomes too large (relative to the typical wave slope) for any of the triads under consideration.

To demonstrate the potential for the simultaneous excitation of two triads within a circular cylinder of finite depth, we consider two scenarios: (i) the excitation of two triads that share a common wave mode; and (ii) the excitation of two triads that do not share any common wave modes. Heuristically, case (ii) is more common than case (i) owing to the number of similar fluid depths in table 1; however, case (i) will likely generate a far richer set of dynamics owing to the nonlinear interaction between the two triads (McEwan, Mander & Smith 1972; Craik 1986; Chow, Henderson & Segur 1996; Choi, Chabane & Taklo 2021). As an example of case (i), we consider triads 11 and 12 in table 1, with nearby critical depths $h_{c,1} = 0.23678$ and $h_{c,2} = 0.21395$, respectively. As mode $(m_3,n_3) = (3,3)$ is common to both triads, inter-triad resonance may arise at an intermediate depth, e.g. $h = 0.225$. Furthermore, an example of case (ii) arises for triads 13 and 14 in table 1, with nearby critical depths $h_{c,1} = 0.19839$ and $h_{c,2} = 0.19814$. Neither of these triads share a common wave mode, so one would not expect the inter-triad energy exchange discussed in case (i). Nevertheless, one might anticipate a signature of these two triads to be visible in the surface evolution for an intermediate depth, e.g. $h = 0.19825$. The theoretical and numerical exploration of coupled triads in a circular cylinder will be the focus of future investigation.

3.3.3. Annular cylinder

A natural variation upon a circular cylinder is an annulus of inner radius $r_0 \in (0,1)$ and outer radius 1. By varying $r_0$, the annulus approaches a circular cylinder as $r_0 \rightarrow 0^+$, and a quasi-one-dimensional periodic ring as $r_0 \rightarrow 1^-$. Notably, resonant triads are impossible for a one-dimensional periodic ring, as can be shown by modifying the arguments presented for the case of a rectangular cylinder (see § 3.3.1). Thus, one might anticipate that the existence of triads in an annular cylinder depends critically on the inner radius, $r_0$. Rather than enumerating some possible triads for given values of $r_0$, we instead track the corresponding critical depth, $h_c$, for the triads identified for a circular cylinder (see table 1) as $r_0$ is progressively increased from zero. Of particular interest is determining whether a given triad exists for all $r_0 < 1$, or whether there is some critical inner radius, $r_c$, beyond which the triad ceases to exist, with either $h_c \rightarrow 0$ or $h_c \rightarrow \infty$ as $r_0 \rightarrow r_c^-$.

The (complex-valued) eigenmodes in an annular domain are cylinder functions of the form

(3.20)\begin{equation} \varPhi_{mn}(r,\theta) = \frac{1}{\mathscr{N}_{mn}}[\mathrm{J}_m(k_{mn}r) \cos(\gamma_{mn}{\rm \pi}) + \mathrm{Y}_m(k_{mn}r) \sin(\gamma_{mn}{\rm \pi})] \mathrm{e}^{\mathrm{i} m \theta}, \end{equation}

where $\mathscr {N}_{mn} > 0$ is a normalisation constant, $\mathrm {Y}_m$ is the Bessel function of the second kind with order $m$ (an integer) and $\gamma _{mn} \in [0,1]$ determines the weighting between the two Bessel functions. As shown in Appendix B, the no-flux condition (see (2.1d)) on the inner and outer walls determines that the wavenumbers, $k_{mn}(r_0)$, satisfy the equation

(3.21)\begin{equation} \mathrm{J}_m'(k_{mn} r_0) \mathrm{Y}_m'(k_{mn}) - \mathrm{J}_m'(k_{mn}) \mathrm{Y}_m'(k_{mn} r_0) = 0. \end{equation}

A formula for the corresponding value of $\gamma _{mn}$ is determined in Appendix B. Once again, the wavenumbers, $k_{mn}$, are ordered so that $0 < k_{m1} < k_{m2} < \cdots$ (excluding $k_{00} = 0$) and satisfy $-\Delta \varPhi _{mn} = k_{mn}^2\varPhi _{mn}$. Three correlated wave modes may form a resonant triad (for a judicious choice of the fluid depth) provided that the corresponding wavenumbers, $K_j$, which depend on the channel width, $1 - r_0$, satisfy the bounds given in Theorem 1.

Bifurcating from the limiting case of a circular cylinder, we track the critical depth (when such a depth exists) of different triads as $r_0$ is progressively increased. The predominant behaviour is characterised by the example presented in figure 3, for which we consider the triad whose critical depth is $h_c = 0.17266$ as $r_0 \rightarrow 0^+$ (see triad 15 in table 1). Given that $h_c$ is fairly small in this limit, one might anticipate that the triad ceases to exist with $h_c \rightarrow 0$; somewhat surprisingly, however, the opposite scenario arises, with $h_c \rightarrow \infty$ as $r_0 \rightarrow r_c^-$ ($r_c \approx 0.57$ in this example). It follows, therefore, that the triad may persist for narrow channels only when the fluid is sufficiently deep. We note, however, that there exist (at least) two relatively rare transitions for increasing $r_0$, which we briefly describe as follows: (i) the triad ceases to exist when $h_c \rightarrow 0$ as $r_0 \rightarrow r_c^-$, which may arise when bifurcating from a sufficiently shallow circular cylinder (e.g. triad 25 in table 1); and (ii) the triad continues to exist for all $r_0 < 1$, with $h_c \rightarrow 0$ and $K_j \rightarrow \infty$ as $r_0 \rightarrow 1$, yet the normalised depth, $h_c K_3$, remains finite, and the normalised wavenumbers, $K_1/K_3$ and $K_2/K_3$, remain within the triad existence region (e.g. triad 14 in table 1). Owing to the appreciable influence of viscous effects for relatively shallow fluids, the physical relevance of these latter two scenarios is somewhat nebulous, however.

Figure 3. The existence and predominant characteristics of a triad in an annular cylinder with inner radius $r_0$ and outer radius 1. The triad bifurcates from the critical depth $h_c = 0.17266$ as $r_0 \rightarrow 0$ (the limiting case of a circular cylinder), with corresponding wavenumbers presented in table 1 (see triad 15). (a) The critical depth, $h_c$ (blue curve), with $h_c \rightarrow \infty$ as $r_0 \rightarrow r_c^-$, where $r_c \approx 0.57$ (black line). (b) The corresponding wavenumbers, $K_j$, all of which remain finite for $r_0 < r_c$ (black line). (c) The normalised wavenumbers, $K_1/K_3$ and $K_2/K_3$, parametrised by increasing $r_0$ (blue arrow), with the limiting case $r_0 \rightarrow 0$ denoted by the white dot. The wavenumbers leave the triad existence region (see Theorem 1) via the left-hand boundary (black curve) as $r_0 \rightarrow r_c^-$.

4. The evolution of resonant triads

Having established the existence of resonant triads, we now determine the long-time triad evolution, utilising the method of multiple scales. Ostensibly, the calculations necessary for determining the triad equations are a variation upon the pioneering work of McGoldrick (1965, 1970a,b) in the absence of surface tension. However, the confinement of the fluid to a cylinder imposes some additional considerations, the salient details of which we outline below. Finally, we note that an alternative approach to multiple scales is Whitham's technique of averaging the system's Lagrangian (Whitham 1965a,b, 1967a,b), which has the advantage of streamlining some algebraic calculations (Simmons 1969; Miles 1976, 1984a); nevertheless, multiple-scale analysis is sufficient for our purposes and allows for the possible inclusion of higher-order corrections in the asymptotic expansion (McGoldrick 1970b).

In a manner similar to § 3, we consider three linear wave modes (with real-valued eigenfunctions), enumerated $n_1$, $n_2$ and $n_3$, where we denote

(4.1a–d)\begin{equation} \varOmega_j = \omega_{n_j}, \quad K_j = k_{n_j}, \quad L_j = \hat{\mathscr{L}}_{n_j}, \quad \mathrm{and}\quad \varPsi_j(\boldsymbol{x}) = \varPhi_{n_j}(\boldsymbol{x}) \quad \mathrm{for}\ j = 1,2,3. \end{equation}

In contrast to § 3, however, we now allow each (non-zero) angular frequency to be either negative or positive: the resonance condition on the angular frequencies is henceforth defined

(4.2)\begin{equation} \varOmega_1 + \varOmega_2 + \varOmega_3 = 0. \end{equation}

The modified requirement on the angular frequencies (4.2) is not restrictive on the possible triad combinations; one may recover (3.4) by mapping $\varOmega _3 \mapsto -\varOmega _3$, for example. The decision behind the summation condition on the angular frequencies is motivated by the cyclical symmetry of (4.2), a property that will be inherited by the resultant amplitude equations (Simmons 1969). As a consequence, one need only derive the amplitude equation for one of the wave modes; the amplitude equations for the remaining two wave modes follow by cyclic permutation of the subscripts $(1,2,3)$.

Before embarking on the multiple-scale analysis presented in § 4.1, we remark upon two caveats. First, we note that (4.2) corresponds to an exact resonance, for which the fluid depth, $h$, is chosen to be precisely equal to the critical depth, $h_c$. In practice, however, there may be a small discrepancy between $h$ and $h_c$, resulting in a the sum of the angular frequencies being slightly offset from zero. When the frequency detuning is sufficiently weak, e.g. $\varOmega _1 + \varOmega _2 + \varOmega _3 = O(\epsilon )$, one may modify the following asymptotic analysis to derive a similar set of amplitude equations (see § 4.3.1). Second, our analysis in § 4.1 is not valid when two of the wave modes coincide. This case corresponds to a 1:2 resonance, for which the corresponding evolution equations were derived by Miles (1976) using Whitham modulation theory (as summarised in § 4.3.2).

4.1. Multiple-scale analysis

In order to determine the evolution of each of the three dominant wave modes involved in an exact resonance, we utilise the method of multiple scales (Kevorkian & Cole 1996; Strogatz 2015). Specifically, we seek a perturbation solution to the Benney–Luke equation (2.3) of the form $u \sim u_0 + \epsilon u_1 + O(\epsilon ^2)$. The leading-order terms in (2.3) determine that $u_0$ satisfies $\partial _{tt}u_0 + \mathscr {L} u_0 = 0$; we choose to consider a leading-order solution comprised only of the three triad modes (all other modes are assumed to be smaller in magnitude and appear at higher order), giving rise to the leading-order form

(4.3)\begin{equation} u_0(\boldsymbol{x} ,t, \tau) = \sum_{j = 1}^3 [ A_j(\tau) \varPsi_j(\boldsymbol{x}) \exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}]. \end{equation}

In (4.3), we have introduced the slow time scale $\tau = \epsilon t$, which governs the evolution of each complex amplitude, $A_j$. As $\epsilon$ and $t$ are both independent variables, we treat $\tau$ and $t$ as independent time scales, giving rise to the transformation of derivatives $\partial _t \mapsto \partial _t + \epsilon \partial _\tau$. Finally, $\mathrm {c.c.}$ denotes the complex conjugate of the preceding term, a contribution necessary for real $u_0$.

So as to determine coupled evolution equations for each complex amplitude, $A_j$, we consider terms of $O(\epsilon )$ in the Benney–Luke equation (2.3). By substituting the leading-order solution, $u_0$, into the nonlinear terms and applying the triad condition for the angular frequencies (4.2), we obtain the following problem for $u_1$:

(4.4)\begin{equation} \partial_{tt}u_1 + \mathscr{L} u_1 ={-}\left[\sum_{j = 1}^3 f_j(\boldsymbol{x},\tau)\exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}\right] + \mbox{non-resonant terms}. \end{equation}

As we will see below, each of the functions $f_j(\boldsymbol {x},\tau )$ appearing on the right-hand side of (4.4) will play a fundamental role when determining the amplitude equations; specifically,

(4.5) \begin{align} f_1 &={-}2\mathrm{i} \varOmega_1 \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} \varPsi_1 + \mathrm{i} A_2^* A_3^*([\varOmega_2 (L_3^2 - K_3^2) + \varOmega_3(L_2^2 - K_2^2) - 2\varOmega_1 L_2L_3]\varPsi_2\varPsi_3\nonumber\\ &\quad - 2\varOmega_1 \boldsymbol{\nabla}\varPsi_2 \boldsymbol{\cdot}\boldsymbol{\nabla} \varPsi_3), \end{align}

where $f_2$ and $f_3$ follow upon cyclic permutation of the subscripts $(1,2,3)$. Finally, we note that the ‘non-resonant terms’ in (4.4) are of the general form $p(\boldsymbol {x},\tau )\mathrm {e}^{\mathrm {i} \varsigma t}$, where we assume that the angular frequency, $\varsigma$, is not equal (or close) to any of the angular frequencies, $\pm \omega _n$, associated with linear wave modes (see § 3).

We proceed by projecting (4.4) onto each of the three wave modes, giving rise to differential equations of the form (for $j = 1,2,3$)

(4.6)\begin{equation} \partial_{tt}\hat{u}_{1,j} + L_j \hat{u}_{1,j} ={-}[\langle \varPsi_j, f_j \rangle \exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}] + \mbox{non-resonant terms}, \end{equation}

where $\hat {u}_{1,j} = \langle \varPsi _j, u_1 \rangle$ is the projection of $u_1$ onto the mode $\varPsi _j$. By recalling that $L_j = \varOmega _j^2$, we immediately see that the term in square brackets in (4.6) is itself a solution to the linear operator $\partial _{tt} + \varOmega _j^2$. It follows that the solution of (4.6) comprises of particular solutions that have temporal dependence $t\exp ({\pm \mathrm {i} \varOmega _j t})$, leading to an ill-posed asymptotic expansion when $\epsilon t = O(1)$. The resolution to this problem is achieved via the solubility condition $\langle \varPsi _j, f_j\rangle = 0,$ which suppresses the secular growth.

By applying the solubility condition $\langle \varPsi _j, f_j\rangle = 0$ for $j = 1,2,3$, we conclude that the complex amplitude, $A_j(\tau )$, of each wave mode, $\varPsi _j(\boldsymbol {x}) \exp ({-\mathrm {i}\varOmega _j t})$, evolves according to the triad system of canonical form (Bretherton 1964; Craik 1986)

(4.7a–c)\begin{equation} \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} = \alpha_1 A_2^* A_3^*, \quad \frac{\mathrm{d} {A_2}}{\mathrm{d} {\tau}} = \alpha_2 A_1^* A_3^*, \quad \frac{\mathrm{d} {A_3}}{\mathrm{d} {\tau}} = \alpha_3 A_1^* A_2^*, \end{equation}

where

(4.8) \begin{equation} \alpha_1 = \frac{1}{2\varOmega_1}( [\varOmega_2 (L_3^2 - K_3^2) + \varOmega_3(L_2^2 - K_2^2) - 2\varOmega_1 L_2L_3] \mathscr{C} - 2\varOmega_1 \langle \varPsi_1, \boldsymbol{\nabla}\varPsi_2 \boldsymbol{\cdot}\boldsymbol{\nabla} \varPsi_3 \rangle), \end{equation}

while $\alpha _2$ and $\alpha _3$ follow by cyclic coefficient of the subscripts $(1,2,3)$. Furthermore, the correlation integral, $\mathscr {C}$, is defined

(4.9)\begin{equation} \mathscr{C} = \frac{1}{S}\iint_\mathcal{D} \varPsi_1\varPsi_2\varPsi_3\,\mathrm{d}A, \end{equation}

where we recall that $S$ is the area of the cylinder cross-section (see § 2.2). As the triad equations (4.7a–c) are valid for $\tau = O(1)$ (or $t = O(1/\epsilon )$), their dynamics yields an informative view of the long-time evolution of the resonant triad.

In order to assess the influence of the triad coefficients on the triad evolution (see § 4.2), we first simplify the algebraic form given in (4.8). As shown by Miles (1976), one may simplify the inner product $\langle \varPsi _1, \boldsymbol {\nabla } \varPsi _2 \boldsymbol {\cdot }\boldsymbol {\nabla } \varPsi _3 \rangle$ by repeated application of the divergence theorem and utilisation of the relationship $-\Delta \varPsi _j = K_j^2\varPsi _j$; it follows that

(4.10)\begin{equation} \langle \varPsi_1, \boldsymbol{\nabla}\varPsi_2 \boldsymbol{\cdot}\boldsymbol{\nabla} \varPsi_3 \rangle = \tfrac{1}{2}(K_2^2 + K_3^2 - K_1^2) \mathscr{C}, \end{equation}

where $\mathscr {C}$ is the correlation integral defined in (4.9). We then substitute (4.10) into (4.8) and simplify using the relation $\varOmega _1 + \varOmega _2 + \varOmega _3 = 0$. After some algebra, we derive the reduced expression

(4.11)\begin{equation} \alpha_1 = \frac{\mathscr{C}}{2\varOmega_1}\left( \varOmega_2 L_3^2 + \varOmega_3L_2^2 - 2\varOmega_1 L_2 L_3 + \sum_{l = 1}^3 \varOmega_l K_l^2 \right), \end{equation}

where $\alpha _2$ and $\alpha _3$ follow similarly. Finally, we demonstrate in Appendix C that the algebraic form of the triad coefficients may be further reduced to

(4.12)\begin{equation} \alpha_j = \frac{\mathscr{C}\beta}{2\varOmega_j} \quad\mathrm{for}\ j = 1,2,3, \end{equation}

where

(4.13)\begin{equation} \beta = \sum_{l= 1}^3 \varOmega_lK_l^2 - \frac{1}{2}\varOmega_1\varOmega_2\varOmega_3 (\varOmega_1^2 + \varOmega_2^2 + \varOmega_3^2). \end{equation}

Equations (4.7a–c), (4.9), (4.12) and (4.13) constitute the triad equations for resonant gravity waves confined to a cylinder of finite depth. Although the triad equations are of canonical form (Bretherton 1964), the novelty of our investigation is the computation of the coefficients, $\alpha _j$, whose algebraic form is specific to our system.

4.2. Properties of the triad coefficients

The simplified form of the coefficients, $\alpha _j$ (4.12), allows for some important theoretical observations that were obfuscated by the more complicated expressions for $\alpha _j$ given in (4.8) and (4.11). In particular, as exactly two of the angular frequencies, $\varOmega _j$, have the same sign, we deduce from (4.12) that the two corresponding coefficients, $\alpha _j$, also have the same sign, with the third coefficient having the opposite sign. By utilising well-known results pertaining to the canonical triad equations, we conclude that all solutions to the triad equations (4.7a–c) are periodic in time, with solutions expressible in terms of elliptic functions (Ball 1964; Bretherton 1964; Simmons 1969; Craik 1986). Typically, these solutions result in an exchange of energy between the comprising modes, although there is a class of periodic solution that, perhaps counter-intuitively, results in zero energy exchange for all time (Case & Chiu 1977; Chabane & Choi 2019). Moreover, it is readily verified that the leading-order energy density, $\mathscr {E}_1 + \mathscr {E}_2 + \mathscr {E}_3$, is conserved, where $\mathscr {E}_j = \varOmega _j^2 |A_j|^2$, consistent with the Hamiltonian structure of the Euler equations (Bretherton 1964; Craik 1986). The reader is directed to the work of Craik (1986) for a more detailed account of the various properties of the canonical triad equations.

Of particular relevance to the evolution of the triad is the quantity $\beta$ (see (4.13)), which, together with $\mathscr {C}$, determines the time scale over which energy exchange arises. In particular, we present the form of $\beta$ in figure 4 for the case $\varOmega _1,\varOmega _2 > 0$ and $\varOmega _3 < 0$. As we will demonstrate below, $\beta < 0$ in this case; in general, the sign of $\beta$ is the same as the sign of the largest (in magnitude) angular frequency, $\varOmega _j$. Notably, $|\beta |$ decreases sharply towards zero as $K_1 + K_2 \rightarrow K_3$, corresponding to the limit $h_c \rightarrow 0$. Similarly, $|\beta |$ approaches zero in the limiting cases $K_1 \ll K_3$ or $K_2 \ll K_3$, corresponding to one low-oscillatory wave mode interacting with two highly oscillatory wave modes. Away from these limiting cases, however, $|\beta |$ depends only weakly on the wavenumbers, $K_j$, suggesting that the correlation integral, $\mathscr {C}$, predominantly controls the time scale of the triad evolution. Finally, we observe that $\beta$ is symmetric about the line $K_1 = K_2$, consistent with the invariance of (4.13) under the mapping $K_1\leftrightarrow K_2$ (and hence, $\varOmega _1 \leftrightarrow \varOmega _2$).

Figure 4. Contours of $\beta/K_3^{5/2}$ (see (4.13)) for the case $\varOmega _1, \varOmega _2 > 0$ and $\varOmega _3 < 0$ (with $\varOmega _1 + \varOmega _2 + \varOmega _3 = 0$), for which $\beta < 0$ (see (4.16)).

We conclude this section by proving that $\beta < 0$ in the case $\varOmega _1, \varOmega _2 > 0$ and $\varOmega _3 < 0$. By comparing the forms of (4.13) and (4.11), and then permuting the subscripts $(1,2,3) \mapsto (3,1,2)$, we first note that $\beta$ may be equivalently expressed as

(4.14)\begin{equation} \beta = \varOmega_1 L_2^2 + \varOmega_2L_1^2 - 2\varOmega_3L_1L_2 + \sum_{l = 1}^3\varOmega_l K_l^2, \end{equation}

or

(4.15)\begin{equation} \beta = \varOmega_1 (L_2^2 + K_1^2) + \varOmega_2(L_1^2 + K_2^2) + |\varOmega_3| ( 2L_1 L_2 - K_3^2). \end{equation}

By bounding $L_j = K_j \tanh (K_j h_c) < K_j$ for $0 < h_c < \infty$ and utilising the relation $\varOmega _1 + \varOmega _2 = |\varOmega _3|$, we obtain

(4.16) \begin{equation} \beta < |\varOmega_3|(K_1^2 + K_2^2 + 2K_1 K_2 - K_3^2) = |\varOmega_3| ((K_1 + K_2)^2 - K_3^2). \end{equation}

As resonant triads exist only when $K_1 + K_2 < K_3$ (see Theorem 1), we conclude that $\beta < 0$ in this case.

4.3. Summary

To summarise our theoretical developments, the velocity potential, $u$, at the fluid rest level ($z = 0$) evolves according to

(4.17)\begin{equation} u(\boldsymbol{x},t) \sim \sum_{j = 1}^3 [A_j(\tau) \varPsi_j(\boldsymbol{x}) \exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}] + O(\epsilon), \end{equation}

where the complex amplitudes, $A_j(\tau )$, evolve over the slow time scale, $\tau = \epsilon t$, according to the triad equations (4.7a–c). In particular, the triad coefficients, $\alpha _j$ (see (4.12)), are defined in terms of the correlation integral, $\mathscr {C}$ (4.9), and the coefficient $\beta$ (4.13). Notably, we assume that $\mathscr {C}$ is non-zero; if this condition were violated then all three of the triad coefficients, $\alpha _j$, would be equal to zero, giving rise to non-interacting wave modes at leading order (contradicting the notion of a triad). Indeed, the condition $\mathscr {C} \neq 0$ is identical to the correlation condition detailed in (3.5), the origins of which we have now justified. Finally, the evolution of the free surface, $\eta$, may be recovered by recalling that $\eta = -u_t + O(\epsilon )$: we conclude that $\eta (\boldsymbol {x},t)$ has a similar leading-order form to $u(\boldsymbol {x},t)$, but each complex amplitude, $A_j(\tau )$, in (4.17) is replaced by $\mathrm {i} \varOmega _j A_j(\tau )$ (see (5.3) below).

We briefly contrast our investigation of triad interaction with the early-time calculation of Michel (2019), who characterised the initial linear growth of a child mode induced by the nonlinear interaction of two parent modes (where all three modes comprise the triad). If modes 1 and 2 are the parent modes and mode 3 is the child mode, then the initial linear growth may be deduced directly from triad equations (4.7a–c) in the limit $|A_3| \ll |A_1| \sim |A_2|$. Specifically, the initial variation of $A_1$ and $A_2$ is slow relative to that of $A_3$, which has the approximate early-time form $A_3(\tau ) \approx \alpha _3 C_1^* C_2^* \tau + C_3$, where $C_j = A_j(0)$. Notably, the linear growth rate of the child mode depends on the corresponding triad coefficient, $\alpha _3$, and the product of the initial amplitudes of the two parent modes. However, our result for circular cylinders differs to that of Michel; we believe that the author neglected some important nonlinear contributions (compare Michel's equation (A2) with (2.4) and (2.4a) of Longuet-Higgins 1962). As Michel's experiment verified the scaling of the interaction only up to a proportionality constant, this discrepancy was not captured.

4.3.1. The influence of weak detuning

As discussed earlier in § 4, the analysis in §§ 4.1 and 4.2 does not account for weak detuning of the angular frequencies, as might arise when the fluid depth, $h$, differs slightly from the critical depth, $h_c$. We now briefly consider the case of weak detuning, for which (4.2) is replaced by the condition $\varOmega _1 + \varOmega _2 + \varOmega _3 = \epsilon \sigma$ (see § 3.3.2); here, $\epsilon$ is the small parameter representative of the typical wave slope (see § 2) and $\sigma = O(1)$ determines the extent of the detuning (Bretherton 1964; McGoldrick 1972). By following a very similar multiple-scale procedure to the case $\sigma = 0$, we obtain amplitude equations that are now augmented by a time-dependent modulation. Specifically, each complex amplitude now evolves according to

(4.18a–c)\begin{equation} \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} = \alpha_1 A_2^* A_3^* \mathrm{e}^{\mathrm{i} \sigma \tau}, \quad \frac{\mathrm{d} {A_2}}{\mathrm{d} {\tau}} = \alpha_2 A_1^* A_3^* \mathrm{e}^{\mathrm{i} \sigma \tau}, \quad \frac{\mathrm{d} {A_3}}{\mathrm{d} {\tau}} = \alpha_3 A_1^* A_2^* \mathrm{e}^{\mathrm{i} \sigma \tau}, \end{equation}

where each coefficient, $\alpha _j$, is defined in (4.12). Although detuning yields non-autonomous amplitude equations, autonomous equations may be derived by mapping $A_j(\tau ) \mapsto A_j(\tau ) \exp ({\mathrm {i} \sigma \tau /3})$ for all $j = 1,2,3$ (Craik 1986). Finally, we note that the energy, $\mathscr {E}_1 + \mathscr {E}_2 + \mathscr {E}_3$, is not exactly conserved when considering the effects of detuning; instead, the energy slowly oscillates about a constant value (Craik 1986).

4.3.2. The case of a 1:2 resonance

A 1:2 resonance is a resonant triad for which two modes comprising the triad coincide. For this case, we define two angular frequencies, $\varOmega _1$ and $\varOmega _2$, so that $\varOmega _2 = 2\varOmega _1$ (Miles 1976), where the connection to resonant triads is clear when writing $\varOmega _1 + \varOmega _1 = \varOmega _2$. By following a very similar multiple-scale procedure to that outlined in § 4.1, we obtain

(4.19)\begin{equation} u(\boldsymbol{x},t) \sim \sum_{j = 1}^2 [A_j(\tau) \varPsi_j(\boldsymbol{x}) \exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}] + O(\epsilon), \end{equation}

where

(4.20a,b)\begin{equation} \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} ={-}\gamma A_1^* A_2\quad \mathrm{and}\quad \frac{\mathrm{d} {A_2}}{\mathrm{d} {\tau}} = \frac{\gamma}{4} A_1^2. \end{equation}

In particular, the evolution of the amplitude equations (4.20a,b) depends on the coefficient $\gamma = \mathscr {C}(K_2^2 - K_1^2 - 3\varOmega _1^4)$, where $\mathscr {C} = ({1}/{S})\iint _{\mathcal {D}} \varPsi _1^2 \varPsi _2 \,\mathrm {d}A$ is the correlation integral. Indeed, the amplitude equations (4.20a,b) and coefficient, $\gamma$, are consistent with the results of Miles (1976) when expressing the evolution of each complex amplitude, $A_j$, in polar form (with appropriate rescaling). Finally, we note that a weak detuning (see § 4.3.1) may also be incorporated within the amplitude equations (4.20a,b), thereby accounting for a slight mismatch between the fluid depth, $h$, and the corresponding critical depth, $h_c$ (Miles 1976).

Of particular interest is the evolution of weakly nonlinear waves steadily propagating around a circular cylinder of unit radius, focusing on the case where the fluid depth is precisely equal to the critical depth of a 1:2 resonance (Yang et al. 2021). For the complex-valued eigenmodes defined in (3.18), the correlation condition, $\iint _{\mathcal {D}} \varPsi _1^2 \varPsi _2^* \,\mathrm {d}A \neq 0$, determines that the angular wavenumbers satisfy $m_2 = 2m_1$ (Chossat & Dias 1995; Yang et al. 2021). By expressing the complex wave amplitudes in polar form, $A_j(\tau ) = a_j(\tau ) \exp ({\mathrm {i} \theta _j(\tau )})$ (for $j = 1,2$), (4.20a,b) may be recast as (Miles 1976)

(4.21a–c)\begin{align} \frac{\mathrm{d} {a_1}}{\mathrm{d} {\tau}} ={-}\gamma a_1 a_2 \cos\varTheta, \quad \frac{\mathrm{d} {a_2}}{\mathrm{d} {\tau}} = \frac{\gamma}{4}a_1^2 \cos\varTheta, \quad \frac{\mathrm{d} {\varTheta}}{\mathrm{d} {\tau}} = 2\gamma a_2 \left[1 - \frac{a_1^2}{8a_2^2}\right] \sin\varTheta, \end{align}

where $\varTheta (\tau ) = \theta _2(\tau ) - 2\theta _1(\tau )$ is the time-dependent phase shift. Steadily propagating waves correspond to time-independent solutions for $a_1$, $a_2$ (both non-zero) and $\varTheta$, from which we deduce that $\cos \varTheta = 0$ and $a_1/a_2 = 2\sqrt {2}$. Indeed, it is remarkable that the amplitude ratio of the two dominant (normalised) wave modes is independent of the angular wavenumbers, $m_j$, the radial wavenumbers, $K_j$, and the corresponding angular frequencies, $\varOmega _j$ (see § 3.3.2 for details). Furthermore, one may readily determine the relationship between the angular velocity of the steady wave rotation and the corresponding wave amplitude, which may then be compared with the numerical solution of the full Euler equations (Yang et al. 2021). This comparison, as well as a comparison with steadily propagating waves computed from various truncations of the Euler equations, will be the subject of future investigation.

5. The excitation of resonant triads

Having established the existence and evolution of resonant triads, we now focus on the excitation of a particular triad via external forcing. So as to motivate the method of excitation, we first recall (§ 5.1) the well-known result that one mode in the triad may, or may not, excite the other two modes (Davis & Acrivos 1967; Hasselmann 1967; Simmons 1969); in the case of excitation, the initial mode is referred to as the pump mode (Craik 1986). We will then utilise the criterion of the pump mode to excite all three modes in the triad via a pulsating pressure source (§ 5.2). Throughout this section, we continue with the convention that the triad angular frequencies satisfy $\varOmega _1 + \varOmega _2 + \varOmega _3 = 0$, as set forth in § 4.

5.1. Excitation via the triad pump mode

To first identify the triad pump mode and then characterise the resultant excitation, we consider the case for which $A_3$, say, is much larger in magnitude than the other two mode amplitudes, so $|A_1|, |A_2| \ll |A_3|$ (Davis & Acrivos 1967; Hasselmann 1967; Simmons 1969). By linearising the triad equations (4.7a–c), we obtain

(5.1a–c)\begin{equation} \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} = \alpha_1 A_2^* A_3^*, \quad \frac{\mathrm{d} {A_2}}{\mathrm{d} {\tau}} = \alpha_2 A_1^* A_3^*, \quad \frac{\mathrm{d} {A_3}}{\mathrm{d} {\tau}} = 0, \end{equation}

from which we immediately conclude that $A_3$ is constant (whilst the linearisation assumption holds); we denote $A_3(\tau ) = C$ for some given complex number $C$. By considering second derivatives of $A_1$ and $A_2$, we deduce the linearised evolution equations (Craik 1986)

(5.2a,b)\begin{equation} \frac{\mathrm{d} {{}^2A_1}}{\mathrm{d} {\tau^2}} = \alpha_1\alpha_2 |C|^2 A_1\quad\mathrm{and}\quad \frac{\mathrm{d} {{}^2A_2}}{\mathrm{d} {\tau^2}} = \alpha_1\alpha_2 |C|^2 A_2, \end{equation}

where $\alpha _1\alpha _2 = \mathscr {C}^2\beta ^2/(4\varOmega _1\varOmega _2)$ (see (4.12)). We conclude that $A_1(\tau )$ and $A_2(\tau )$ grow exponentially in time (whilst the linearisation approximation holds) when $\varOmega _1 \varOmega _2 > 0$, and exhibit sinusoidal oscillations when $\varOmega _1\varOmega _2 < 0$ (Davis & Acrivos 1967; Hasselmann 1967; Craik 1986). Thus, mode 3 may excite modes 1 and 2 when $\varOmega _1$ and $\varOmega _2$ have the same sign (and likewise for other mode permutations). As one angular frequency must have a different sign from the other two (so as to satisfy $\varOmega _1 + \varOmega _2 + \varOmega _3 = 0$), we conclude that the mode whose angular frequency is largest in magnitude (i.e. differs in sign) is the triad pump mode (Craik 1986). Equivalently, the pump mode is the mode with largest wavenumber, $K_j$, providing a robust mechanism for an inverse energy cascade to lower wavenumbers (Annenkov & Shrira 2006).

To visualise the influence of the pump mode on the resultant free-surface pattern, we present the solution of the triad equations (4.7a–c) and the corresponding pump-mode approximation (5.1a–c) in figure 5. By recalling that the free surface satisfies $\eta = -u_t + O(\epsilon )$, we first deduce that

(5.3)\begin{equation} \eta(\boldsymbol{x},t) \sim \sum_{j = 1}^3 [\mathrm{i} \varOmega_j A_j(\tau) \varPsi_j(\boldsymbol{x}) \exp({-\mathrm{i} \varOmega_j t}) + \mathrm{c.c.}] + O(\epsilon). \end{equation}

For the case of a circular cylinder, we utilise the complex-valued eigenmodes defined in (3.18), corresponding to the superposition of steadily propagating waves for $m_j \neq 0$ (the rotation direction depends on the sign of $\varOmega _j /m_j$). Upon initialising the system so that the energy is primarily within the pump mode (mode 3), modes 1 and 2 are gradually excited due to nonlinear interaction, with exponential growth evident for $\tau \lesssim 10$. As time further increases, the dynamics departs from the pump-mode approximation: the energy in the pump mode appreciably decreases, whilst the energy in modes 1 and 2 saturates. The free surface varies qualitatively during this evolution, with an appreciable change in pattern structure visible by $\tau = 24$ (primarily a superposition of modes 1 and 2). Notably, the system evolution is periodic, which becomes apparent over longer time scales.

Figure 5. Excitation of a triad via its pump mode for the case of a circular cylinder. We consider triad 11 in table 1, but with $m_3 \mapsto -m_3$. We choose $\varOmega _1, \varOmega _2 > 0$ and $\varOmega _3 < 0$, so that mode 3 is the pump mode. (a) Evolution of the free surface, $\eta \sim -u_t$, over the slow time scale, $\tau = \epsilon t$, with $\epsilon = 10^{-3}$. (b) The evolution of the wave amplitudes, $|A_j|$, according to the triad equations ((4.7a–c), solid curves) and the pump-mode approximation ((5.1a–c), dashed-dotted curves). Insets: modes 1 (blue), 2 (red) and 3 (gold) at $\tau = 0$; all three modes rotate counter-clockwise. The simulations were initialised from $A_1(0) = 0.01$ and $A_2(0) = 0.01\mathrm {i}$, where $A_3(0)$ was chosen to be the positive real number satisfying $\mathscr {E}_1 + \mathscr {E}_2 + \mathscr {E}_3 = 1$, with $\mathscr {E}_j = \varOmega _j^2 |A_j|^2$ (see § 4.2).

5.2. Excitation via an applied pressure source

Based on the ideas of the previous section, we consider a methodology for exciting the pump mode of a triad, which will subsequently excite the remaining two modes (provided that the initial disturbance of each of the remaining modes is non-zero). Notably, several methods for exciting internal resonances have been considered in prior investigations, primarily focusing on imposed motion of the fluid vessel via horizontal (Miles 1976, 1984c) or vertical vibration (Miles 1976, 1984b; Miles & Henderson 1990; Henderson & Miles 1991). Furthermore, one may, in principle, utilise sinusoidal paddles or plungers to excite a particular triad's pump mode for a given geometry (similar wave makers are used in rectangular wave tanks McGoldrick 1970a; Henderson & Hammack 1987). However, for large-scale fluid tanks, imposed motion of the vessel may be impractical (if the tank were set in a concrete base, for example), and it may be challenging to determine the correct paddle motion necessary to excite a chosen pump mode for geometrically complex cylinders. We choose, therefore, to consider a slightly different approach: we instead excite the pump mode via a pulsating pressure source located just above the free surface (e.g. an air blower).

In order to incorporate a pressure source within our mathematical framework, we first reformulate the dimensionless dynamic boundary condition (2.1b) as

(5.4)\begin{equation} \phi_t + \eta + \frac{\epsilon}{2}(|\boldsymbol{\nabla} \phi|^2 + \phi_z^2) + \epsilon P(\boldsymbol{x}, t) = 0 \quad \mbox{for}\ \boldsymbol{x} \in \mathcal{D}, \quad z = \epsilon \eta, \end{equation}

where the dimensional pressure is $\epsilon ^2 a \rho g P$ for fluid density $\rho$ ($P = 0$ corresponds to atmospheric pressure). The pressure source is chosen to be small in magnitude so that the resultant wave excitation arises over the slow time scale, $\tau = \epsilon t$, and may thus be saturated by weakly nonlinear effects. By modifying the developments outlined in § 2.1, we derive the forced Benney–Luke equation

(5.5)\begin{equation} u_{tt} + \mathscr{L} u + \epsilon\left( u_t(\mathscr{L}^2 + \varDelta)u + \frac{\partial {}}{\partial {t}}[(\mathscr{L} u)^2 + |\boldsymbol{\nabla} u|^2] + \partial_t P\right) = O(\epsilon^2) \quad\mathrm{for}\ \boldsymbol{x} \in \mathcal{D}, \end{equation}

which will be the starting point for the asymptotic analysis.

Before proceeding further, we first describe two forms of the pressure source relevant to our investigation. For a stationary pressure source oscillating periodically over the fast time scale, $t$, we express $P(\boldsymbol {x},t) = f(\tau ) s(\boldsymbol {x}) \exp ({-\mathrm {i} \varOmega _p t}) + \mathrm {c.c.}$, where $s(\boldsymbol {x})$ is a fixed spatial profile (generally spanning the cavity), $f(\tau )$ accounts for a slow modulation in the magnitude of the pressure, and $\varOmega _p$ is the pulsation angular frequency. We choose $\varOmega _p$ to be close to the angular frequency of the pump mode, which, without loss of generality, we assume to be mode 3 (i.e. $\varOmega _3$ has the opposite sign from $\varOmega _1$ and $\varOmega _2$). We denote, therefore, $\varOmega _p = \varOmega _3 + \epsilon \mu$, where $\mu = O(1)$ determines the extent of the frequency mismatch. For a pressure source orbiting the centre of a circular cylinder at a constant angular velocity, we instead posit that $P$ has the form $P(r,\theta, t) = f(\tau ) s(r, \theta - \varOmega _p t)$, where $\varOmega _p = (\varOmega _3 + \epsilon \mu )/m_3$ is the angular velocity of the pressure source (assuming that the pump mode is non-axisymmetric, i.e. $m_3 \neq 0$).

For both standing and orbiting pressure sources, we now follow a similar multiple-scale procedure to that outlined in § 4.1, starting from the forced Benney–Luke equation (5.5). So as to discount the possibility that the pressure source excites more than one mode in the triad, we assume that neither $|\varOmega _1|$ or $|\varOmega _2|$ are close to $|\varOmega _3|$. Furthermore, we incorporate a weak detuning in the triad angular frequencies, denoting $\varOmega _1 + \varOmega _2 + \varOmega _3 = \epsilon \sigma$ (see § 4.3.1). It follows that each complex amplitude, $A_j(\tau )$, evolves according to

(5.6a–c)\begin{equation} \left.\begin{gathered} \frac{\mathrm{d} {A_1}}{\mathrm{d} {\tau}} = \alpha_1 A_2^* A_3^* \mathrm{e}^{\mathrm{i} \sigma \tau}, \quad \frac{\mathrm{d} {A_2}}{\mathrm{d} {\tau}} = \alpha_2 A_1^* A_3^* \mathrm{e}^{\mathrm{i} \sigma \tau}, \\ \frac{\mathrm{d} {A_3}}{\mathrm{d} {\tau}} = \alpha_3 A_1^* A_2^* \mathrm{e}^{\mathrm{i} \sigma \tau} - \varOmega_3s_3f(\tau)\exp({-\mathrm{i} \mu \tau}), \end{gathered}\right\} \end{equation}

where the coefficients, $\alpha _j$, are defined in (4.12). Notably, the pump mode may only be excited provided that the corresponding eigenmode is non-orthogonal to the pressure source, corresponding to a non-zero projection, i.e. $s_3\neq 0$, where $s_3 = \langle \varPsi _3, s\rangle$. Similar equations describing the evolution of forced resonant triads have been explored by McEwan et al. (1972) (with the inclusion of linear damping) and Raupp & Silva Dias (2009).

In the special case of time-independent forcing ($\,f$ constant) and no frequency detuning ($\sigma = \mu = 0$), the dynamics of the forced triad equations has been analysed by Harris, Bustamante & Connaughton (2012), with both periodic and quasi-periodic dynamics reported. We also consider this case, leaving the effects of detuning and variable forcing for future investigation. In this setting, when $|A_1|$, $|A_2|$ and $|A_3|$ are initially small relative to the magnitude of the forcing, $|\varOmega _3 s_3 f|$, the initial growth in $A_3$ is approximately linear (see figure 6a). As mode 3 is the pump mode, the growth in $A_3$ excites $A_1$ and $A_2$, thus activating the triad. The conservation laws of the forced triad equations (Harris et al. 2012) result in a temporary diminution of mode 3, which is later augmented by the external forcing; whence the process repeats. In some parameter regimes, the resulting evolution of the forced triad is periodic in time (see figure 6(b) and Raupp & Silva Dias 2009); in contrast to the findings of Harris et al. (2012), however, we also identify initial conditions (with all other parameters unchanged) that result in hitherto unidentified chaotic dynamics (see figure 6c).

Figure 6. Evolution of the forced triad equations (5.6a–c) for $\sigma = \mu = 0$ and constant $f$. We consider the same triad as figure 5, with $s_3 f = 0.1$. In all four panels, $A_1(0) = 0.02\mathrm {i}$ and $A_2(0) = 0.01$. For $A_3(0) = 0.01$, we observe (a) the initial excitation of the triad and (b) the resultant periodic dynamics (the initial growth is highlighted within the grey box). (c) For $A_3(0) = 0.01\mathrm {i}$, the triad evolution appears chaotic. (d) The separation distance, $\delta (\tau )$, of a trajectory randomly perturbed at time $\tau _{{pert}}$, with $\delta (\tau _{{pert}}) = 10^{-10}$ and the same initialisation (at $\tau = 0$) as (c). The black curves bound $\delta (\tau )$ for 30 equally spaced values of $\tau _{{pert}}$ in the interval $100 \leq \tau _{{pert}} \leq 115$; the blue curve corresponds to $\tau _{{pert}} = 100$.

To verify the chaotic nature of this latter example, we consider the separation distance of two initially adjacent trajectories in phase space, for which we observe exponential divergence in time (see figure 6d). This exponential divergence is indicative of a positive maximal Lyapunov exponent (Strogatz 2015), which characterises the sensitivity to initial conditions exhibited by chaotic systems. Specifically, for a solution, $A_j(\tau )$, and its perturbation, $A_{j,{pert}}(\tau )$, we consider the evolution of the separation distance, defined as

(5.7)\begin{equation} \delta(\tau) = \sqrt{\sum_{j = 1}^3 |A_j(\tau) - A_{j,{pert}}(\tau)|^2 }. \end{equation}

To compute $\delta (\tau )$, we first simulate the forced triad equations (5.6a–c) on the interval $0 \leq \tau \leq \tau _{{pert}}$, with $A_j(0) = A_{j,{pert}}(0)$ and $\tau _{{pert}}$ chosen to be sufficiently large so as to ensure that the chaotic attractor (should one exist) be approached. At $\tau = \tau _{{pert}}$, the real and imaginary parts of $A_{j,{pert}}$ are both randomly perturbed according to a uniform distribution on the interval $(-1, 1)$, with the resulting complex perturbation scaled so that $\delta (\tau _{{pert}}) = 10^{-10}$. We then evolve $A_j(\tau )$ and $A_{j,{pert}}(\tau )$ for $\tau \geq \tau _{{pert}}$, giving rise to initial exponential growth of $\delta (\tau )$, with saturation when the perturbation distance is comparable to the ‘diameter’ of the chaotic attractor (see figure 6d).

To confirm that the exponential growth was not specific to a perturbation about a particular point on the chaotic attractor (Strogatz 2015), we considered 30 equally spaced values of $\tau _{{pert}}$ on the interval $100 \leq \tau _{{pert}} \leq 115$, roughly corresponding to the time taken for one ‘loop’ of the chaotic attractor to take place (see figure 6c). Each simulation was computed with a fourth-order Runge–Kutta method and a time step of 0.005. For each value of $\tau _{{pert}}$, we observed similar exponential divergence of trajectories (see figure 6d); moreover, the evolution of $\delta (t)$ during the growth phase remained unchanged when the time step was decreased to 0.001, with numerical errors only accumulating over longer time scales. Our results thus provide strong evidence that there is a positive maximal Lyapunov exponent in this particular portion of parameter space, indicative of a chaotic dynamics.

6. Discussion

We have performed a systematic investigation into nonlinear resonant triads of free-surface gravity waves confined to a cylinder of finite depth; previously studied 1:2 resonances are obtained as special cases. A key result of our study is Theorem 1, which determines whether there exists a fluid depth at which three given wave modes resonate due to the nonlinear evolution of the fluid. Equipped with this result, we determined the long-time fluid evolution using multiple-scale analysis, from which we deduced that all solutions to the triad equations are periodic in time. Finally, we determined that a given triad may be excited via external forcing of the triad's pump mode, thereby providing a mechanism for exciting a given triad in a wave tank. All our results are derived for cylinders of arbitrary cross-section (barring some technical assumptions; see § 2), thus forming a broad framework for characterising nonlinear resonance of confined free-surface gravity waves. In particular, our theoretical developments buttress experimental observations (Michel 2019) and demonstrate the potential generality of confinement as a mechanism for promoting nonlinear resonance.

A second fundamental component of our study is the influence of the cylinder cross-section on the existence of resonant triads; for example, resonant triads are impossible in rectangular cylinders, yet abundant within circular and annular cylinders (for particular fluid depths). Of the vast array of resonances arising in a circular cylinder, triads consisting of an axisymmetric pump mode and two identical counter-propagating waves are of notable interest. This combination of axisymmetric and non-axisymmetric modes possesses an interesting analogy to the excitation of counter-propagating subharmonic beach edge waves due to a normally incident standing wave (Guza & Davis 1974). Specifically, the wave crests of the standing axisymmetric mode are always parallel to the bounding wall of the circular cylinder, and may excite steadily propagating waves that are periodic in the azimuthal direction. For the special case for which the amplitudes of the two counter-propagating modes coincide, one observes the resonant interaction of standing axisymmetric and non-axisymmetric waves.

So as to gain a deeper insight into the influence of nonlinearity on resonant triads, a primary focus for future investigations will be the simulation of the Euler equations within a cylindrical domain, with consideration of various truncated systems (Craig & Sulem 1993; Milewski & Keller 1996; Berger & Milewski 2003; Wang & Milewski 2012). From a computational perspective, the most natural geometry to consider is a circular cylinder (Qadeer & Wilkening 2019); this geometry has been previously explored in the context of steadily propagating nonlinear waves in the vicinity of a 1:2 resonance (Bryant 1989; Yang et al. 2021), but it remains to assess the efficacy of the amplitude equations (4.7a–c) for predicting the evolution of nonlinear triads. Indeed, exploration of the nonlinear dynamics may reveal additional resonant triads arising beyond the small-wave-amplitude limit explored herein. Of similar interest is the fluid evolution when multiple triads are excited at a single depth, with the potential for energy exchange via triad–triad interactions (McEwan et al. 1972; Craik 1986; Chow et al. 1996; Choi et al. 2021). The simulation of free-surface gravity waves in non-circular cylinders presents a more formidable challenge, however, except for cylinder cross-sections that possess a tractable eigenmode decomposition.

A second natural avenue for future investigation is to characterise the influence of applied forcing on resonant triads. For example, when the fluid bath is subjected to sufficiently vigorous vertical vibration, Faraday waves (Faraday 1831; Kumar 1996) may appear on the free surface; although this scenario has been studied in the case of a 1:2 internal resonance (Miles 1984b; Miles & Henderson 1990; Henderson & Miles 1991), resonant triads may give rise to the formation of more exotic free-surface patterns, particularly at fluid depths that differ from that of a 1:2 resonance. In a similar vein, horizontal vibration (Miles 1976, 1984c) or a pulsating pressure source at the frequency of the triad's pump mode may lead to a wealth of periodic and quasi-periodic dynamics, as predicted by the forced triad equations (Harris et al. 2012). Our study has indicated, however, that a chaotic dynamics is also possible in some parameter regimes, and might thus be excited in numerical simulation or experiments. Lastly, our study has focused on flat-bottomed cylinders; it seems plausible, however, that submerged topography may enhance or mitigate certain resonances, which may be an important consideration in the design of industrial-scale fluid tanks.

Finally, our study has focused on the special case of a liquid–air interface, for which the dynamics of the air is neglected within the Euler equations. It is natural, however, to extend our formulation to the case of two-layer flows (in the absence of surface tension), with two immiscible fluids (e.g. air and water) confined within a cylinder whose lid and base are both rigid. In this setting, the density difference across the fluid–fluid interface has a strong influence of the system dynamics; it seems plausible, therefore, that additional resonances may be excited in this configuration, relative to the liquid–air interface considered herein. Notably, the anticipated resonances would arise across a single interface, in contrast to the cross-interface resonances explored in previous investigations (Ball 1964; Simmons 1969; Joyce 1974; Segur 1980; Taklo & Choi 2020; Choi et al. 2021). Finally, exploring the influence of parametric forcing (Kumar & Tuckerman 1994) on resonant triads arising for two-layer flows opens up exciting new vistas in nonlinear resonance induced by confinement.