1. Plane water waves
To fix the ideas, we will focus on inviscid gravity waves in the deep-water setting known to lead to a nonlinear Schrödinger equation (NLS) in the approximation of a narrow wavepacket centred around a wavenumber 
$k_{0}$. The latter together with the wave amplitude 
$a$ and gravitational acceleration g set the natural scales for coordinates 
$\boldsymbol {x} \rightarrow k_{0}^{-1} \boldsymbol {x}$, time 
$t \rightarrow \omega _{0}^{-1} t$, surface deflection from the quiescent state 
$\eta \rightarrow a \eta$, and velocity potential 
$\phi \rightarrow (a \omega _{0} / k_{0}) \phi$, where the circular frequency 
$\omega _{0} = \omega (k_{0})$ is dictated by the dispersion relation 
$\omega _{0}^{2} = g k_{0}$ and the factor 
$a \omega _{0} / k_{0}$ follows from balancing the fluid velocity at the free surface with that of the surface itself, 
$\phi _{z} \sim \eta _{t}$. The resulting potential-flow system reads
$$\begin{gather} z \leqslant \varepsilon \eta : \quad \Delta \phi = 0 \quad \text{with} \ \boldsymbol{\nabla}\phi \rightarrow 0, \ z \rightarrow - \infty, \end{gather}$$
$$\begin{gather}z = \varepsilon \eta : \quad \phi_{z} = \eta_{t} + \varepsilon \boldsymbol{\nabla}_{{\perp}} \phi \boldsymbol{\cdot} \boldsymbol{\nabla}_{{\perp}} \eta \quad \text{and} \quad \phi_{t} + \eta + \frac{\varepsilon}{2} |\boldsymbol{\nabla} \phi|^{2} = 0, \end{gather}$$
where the operators will be specialized either to plane 
$\boldsymbol {x}=(x,z)$ in the present section or cylindrical 
$\boldsymbol {x}=(r,z)$ coordinates in § 2 with 
$\boldsymbol {\nabla } = (\boldsymbol {\nabla }_{\perp },\partial _{z})$. The sole parameter 
$\varepsilon = a k_{0}$ in the kinematic (1.1b) and dynamic (1.1c) boundary conditions (BCs) is treated as small since we are interested in reduced-order models compared with (1.1). Linearization of the latter, following the Stokes (Reference Stokes1847) idea, produces
$$\begin{gather} z \leqslant 0: \quad \Delta \phi = 0 \quad \text{with} \ |\boldsymbol{\nabla} \phi| \rightarrow 0, \ z \rightarrow - \infty, \end{gather}$$
$$\begin{gather}z = 0: \quad \phi_{z} + \varepsilon \phi_{zz} \eta + \frac{\varepsilon^{2}}{2} \phi_{zzz} \eta^{2} = \eta_{t} + \varepsilon \boldsymbol{\nabla}_{{\perp}} \eta \boldsymbol{\cdot} \left(\boldsymbol{\nabla}_{{\perp}} \phi + \varepsilon \boldsymbol{\nabla}_{{\perp}} \phi_{z} \eta\right), \end{gather}$$
$$\begin{gather}\phi_{t} +\varepsilon \phi_{tz} \eta + \frac{\varepsilon^{2}}{2} \phi_{tzz} \eta^{2} + \eta ={-} \frac{\varepsilon}{2} \left[(\boldsymbol{\nabla}_{{\perp}} \phi)^{2} + \phi_{z}^{2} + 2 \varepsilon \left(\boldsymbol{\nabla}_{{\perp}} \phi \boldsymbol{\cdot} \nabla_{{\perp}} \phi_{z} + \phi_{z} \phi_{zz}\right) \eta\right], \end{gather}$$
where we kept terms up to 
${O}(\varepsilon ^{2})$. We will look for solutions in the series 
$\phi = \phi _{0} + \varepsilon \phi _{1} + \varepsilon ^{2} \phi _{2} + \ldots$, 
$\eta = \eta _{0} + \varepsilon \eta _{1} + \varepsilon ^{2} \eta _{2} + \ldots$, at each order 
$i$ requiring 
$|\boldsymbol {\nabla } \phi _{i}| \rightarrow 0$, 
$z \rightarrow - \infty$; for the boundary values we will use the same notation, e.g. the free-surface potential 
$\phi (t,x,0) \equiv \phi (t,x)$.
1.1. Standard multiple-scale view
Solving (1.1) with such regular perturbations proves to lead to secular divergencies (Hakim Reference Hakim1998), which necessitates the introduction of multiple scales 
$(t,x,z) \rightarrow (t,T = \varepsilon t, \tau = \varepsilon ^{2} t; x, X = \varepsilon x; z, Z = \varepsilon z)$, the origin of which can be gleaned from an 
${O}(\varepsilon ^{0})$-solution of (1.2a) for the right-travelling wave
\begin{equation} \phi_{0}(t,x,z) = \frac{1}{2{\rm \pi}} \int_{0}^{+\infty}{\hat{\phi}_{0}(k) \exp({|k| z})\exp[{\mathrm{i} (k x - \omega t)}] \,\mathrm{d} k} + \mathrm{c.c.}, \quad \omega(k) = |k|^{1/2}, \end{equation}
with the initial boundary data 
$\hat {\phi }_{0}(k)$ being the Fourier transform of 
$\phi _{0}(0,x,0)$ related to the initial free-surface deflection 
$\eta _{0}(0,x)$ via 
$\hat {\eta }_{0}(k) = \mathrm {i} \omega (k) \hat {\phi }_{0}(k)$. Starting the evolution with a narrow wavepacket 
$k=1+\varepsilon \kappa$, (1.3) can be expanded to yield 
$\phi _{0}(t,x,z) \sim \varepsilon \varPhi _{0}(T,\tau,X,Z) \exp ({\mathrm {i} (x - t) + z}) + \mathrm {c.c.}$ with the amplitude
\begin{equation} \varPhi_{0}(T,\tau,X,Z) = \frac{1}{2{\rm \pi}} \int_{-\infty}^{+\infty}{\hat{\phi}_{0}(\kappa) \exp\left({\mathrm{i} \left[\kappa (X - \omega^{\prime}_{0} T) - \frac{1}{2} \omega^{\prime\prime}_{0} \kappa^{2} \tau\right] + \kappa Z}\right) \mathrm{d} \kappa}, \end{equation}
where 
$\omega _{0} = \omega (1) = 1$, 
$\omega ^{\prime }_{0} = \frac {1}{2}$ and 
$\omega ^{\prime \prime }_{0} = - \frac {1}{4}$. Differentiating (1.4) produces
\begin{equation} \varPhi_{0T} + \omega^{\prime}_{0} \varPhi_{0 X} = 0, \quad \mathrm{i} \varPhi_{0 \tau} + \tfrac{1}{2} \omega^{\prime\prime}_{0} \varPhi_{0 XX} = 0 \quad \text{at} \ Z=0, \end{equation}
with the former stating that the wavepacket propagates at the group velocity 
$\omega _{0}^{\prime }$ and the latter being the linear part of the NLS in the laboratory frame of reference. For the subsequent discussion, we briefly review the multiple-scale derivation of the NLS also highlighting how (1.5) formally emerge. At the leading order, the system
\begin{equation} {O}(\epsilon^{0}): \left\{\phi_{0 zz} + \phi_{0 xx} = 0 \ z < 0\ | \phi_{0 z} - \eta_{0 t} = 0 \quad \text{and} \quad \phi_{0 t} + \eta_{0} = 0 \quad \text{at} \ z = 0 \right. \end{equation}
admits solution (1.3) for any 
$k$, in particular, a right-travelling narrow wavepacket 
$\phi _{0} = \varPhi _{0}(T,\tau,X,Z) \exp ({\mathrm {i} (x - t) + z}) + \mathrm {c.c.}$. At the next order, the entire right-hand side of the Poisson equation
\begin{equation} {O}(\epsilon^{1}): \quad \phi_{1 zz} + \phi_{1 xx} ={-} 2 \left(\phi_{0 zZ} + \phi_{0 xX}\right) \quad \text{for} \ z \leqslant 0 \end{equation}
produces secular terms containing exponents 
$\exp [{\pm \mathrm {i} (x - t)}]$, the factors of which vanish provided that
\begin{equation} \varPhi_{0 Z} + \mathrm{i} \varPhi_{0 X} = 0 \quad \text{for} \ Z \leqslant 0, \end{equation}
that is, 
$\varPhi _{0}$ propagates along the characteristic of the first-order elliptic operator 
$\mathrm {i} \partial _{X} + \partial _{Z}$. Although the slow 
$Z$-dependence has been the subject of some debate (Mei, Stiassnie & Yue Reference Mei, Stiassnie and Yue2005), its cause is transparent from (1.3) and disregarding it in (1.8) would otherwise lead to spurious transient growth of the 
$\phi _{1}$-solution in the vertical direction. The combined BC, derived by eliminating 
$\eta _{1}$ from the corresponding kinematic and dynamic BCs, is
\begin{align} {O}(\epsilon^{1}):\quad \phi_{1 z} + \phi_{1 tt} = &- 2 \phi_{0 t T} - \phi_{0 Z} - \phi_{0 z t} \left(\eta_{0 t} + \phi_{0 z}\right) \nonumber\\ &- \eta_{0} \left(\phi_{0 z tt} + \phi_{0 zz}\right) + \phi_{0 x} \left(\eta_{0 x} - \phi_{0 x t}\right) \quad \text{at} \ z = 0. \end{align}
Even though we started with initial conditions (ICs) being a wavepacket at 
$k=1$ and zero elsewhere, nonlinear interactions on the right-hand side of (1.9) generate a clustered mode-distribution 
$\exp [{\mathrm {i} n (x - t)}]$, 
$n \in \Bbb {Z}$ (Eckhaus Reference Eckhaus1993). However, substitution of (1.8) evaluated at 
$Z=0$ and the 
${O}(\epsilon ^{0})$-solutions brings out the solvability condition (1.5a) thus nullifying the linear terms on the right-hand side of (1.9) and enforcing 
$\omega ^{\prime }_{0}(1) = \frac {1}{2}$, whereas the nonlinear terms vanish identically and hence the multiple harmonics are not excited in 
$\phi _{1}$. This results in the main harmonic 
$\phi _{1} = \varPhi _{1}(T,\tau,X,Z) \exp ({\mathrm {i} (x - t) + z}) + \mathrm {c.c.}$, which is identical to that for 
$\phi _{0}$ and therefore can be put to zero, while 
$\eta _{1}$ proves to be non-trivial and, as follows from the dynamic BC at this order, contains the main and double modes.
At the second order, the Poisson equation
\begin{equation} {O}(\epsilon^{2}): \quad \phi_{2 zz} + \phi_{2 xx} ={-} \left(\phi_{0 ZZ} + \phi_{0 XX}\right) \quad \text{for} \ z \leqslant 0 \end{equation}
does not contain any secular terms because the right-hand side vanishes in view of (1.8), i.e. 
$\varPhi _{0}$ is an analytical function in the 
$(X,Z)$-plane, 
$\varPhi _{0 ZZ} + \varPhi _{0 XX}=0$. The combined kinematic and dynamic BCs (not shown due to excessive number of terms) at this order, after simplifying with (1.5), (1.8) and transforming to the 
$\xi =X - \omega _{0}^{\prime } T$ frame, lead to the solvability condition in the form of the NLS
\begin{equation} \mathrm{i} \varPhi_{0 \tau} - \tfrac{1}{8} \varPhi_{0 \xi\xi} - 2 |\varPhi_{0}|^{2} \varPhi_{0} = 0 \quad \text{at} \ Z=0, \end{equation}
i.e. (1.5b) is the linear part of the NLS producing yet another Taylor expansion coefficient 
$\omega ^{\prime \prime }_{0}(1) = - \frac {1}{4}$. Notably, 
$\eta _{1}$ containing main and double harmonics does not contribute to the envelope equation (1.11).
1.2. Hilbert wave equation
In the foregoing multiple-scale approach, we had to solve the problem at each order both in the bulk and at the free surface, but our ultimate focus is the dynamics at the free surface only, i.e. the envelope equation (1.11). Also, the above standard analysis does not reveal the true nature of the ‘wave’ equation for the free-surface dynamics. To comprehend the latter, we may appeal to the solution 
$\phi (t,x,z) = - ({1}/{{\rm \pi} })\int _{-\infty }^{+\infty }{\phi _{x^{\prime }}(t,x^{\prime }) \operatorname {atan}[({x-x^{\prime }})/{z}] \, \mathrm {d} x^{\prime }}$ of the Laplace equation (1.2a) with the BC 
$\phi _{x}(t,x,0)=\phi _{x}(t,x)$ enabled by the fact that the domain for (1.2a) is the lower half-plane and hence constructed from the Poisson-Schwarz integral formula. Differentiating 
$\phi$ leads to a representation in terms of the Hilbert transform (HT) 
$\mathcal{H}$:
\begin{equation} \phi_{z}(t,x,0) = \frac{1}{\rm \pi} {\int\hskip -1,05em -\,}_{-\infty}^{+\infty}{\frac{\phi_{x^{\prime}}(t,x^{\prime})}{x-x^{\prime}} \, \mathrm{d} x^{\prime}} \equiv \mathcal{H}\phi_{x}(t,x) = K * \phi_{x} ={-} \frac{1}{\rm \pi} {\int\hskip -1,05em -\,}_{-\infty}^{+\infty}{\frac{\phi(t,x^{\prime})}{(x-x^{\prime})^{2}} \, \mathrm{d} x^{\prime}}; \end{equation}
here 
${\int\hskip -1,05em -\,}$ stands for the Cauchy principal value integral and by 
$K(x) = {1}/{{\rm \pi} x}$ we denote the Cauchy kernel. While naming the transform after Hilbert by G.H. Hardy has been questioned, cf. discussion in King (Reference King2009), 
$\phi (t,x,z)$ represents a solution of the Riemann–Hilbert problem (Muskhelishvili Reference Muskhelishvili1953) of finding a function 
$\chi (t,x,z) = \phi + \mathrm {i} \psi$ holomorphic in the lower half-plane from the known boundary value 
$\chi (t,x)$ at 
$z=0$. The latter satisfies 
$\mathcal {H}\chi = \mathrm {i} \chi$ per the Sokhotskii–Plemelj theorem, whence it is straightforward to deduce (1.12) using the Cauchy–Riemann equations, thus once again justifying the term coined by Hardy, as the Hilbert's contribution is at the heart of the story, and illuminating the meaning of HT.
Even though holomorphic functions may, in general, behave wildly near boundaries (Mashreghi Reference Mashreghi2022), function (1.3) is analytic up to and including the free surface, which enables linearization (1.2) and allows us to find higher-order derivatives normal to the free surface, needed in (1.2b,c), directly from the Laplace equation (1.2a). Eliminating 
$\eta (t,x)$ from (1.2b,c) and invoking the property 
$\mathcal {H}\phi _{x} = \mathcal {H}_{x}\phi$ furnishes an integro-differential equation for 
$\phi (t,x)$:
\begin{align} \phi_{tt} +
\mathcal{H}_{x}\phi &=\varepsilon
\left\{\hskip1pt\partial_{t}\left[\phi_{t} \mathcal{H}_{tx}\phi -
\frac{1}{2} \left(\phi_{x}^{2} +
(\mathcal{H}_{x}\phi)^{2}\right)\right] -
\partial_{x}\left(\phi_{t}\phi_{x}\right)\right\}
\nonumber\\
&\quad\varepsilon^{2}\left\{\partial_{t}\left[\phi_{t}
\phi_{x}^{2}
\partial_{x}\!\left(\frac{\mathcal{H}_{x}\phi}{\phi_{x}}\right)
+ \frac{1}{2} \phi_{t}^{2} \phi_{txx} - \phi_{t}
(\mathcal{H}_{tx}\phi)^{2} - \frac{1}{2}\!
\left(\phi_{x}^{2} + (\mathcal{H}_{x}\phi)^{2}\right)
\mathcal{H}_{tx}\phi\!\right]\right. \nonumber\\
&\quad\quad \ \ +\phi_{xx} \left[\phi_{t} \mathcal{H}_{tx}\phi -
\frac{1}{2} \left(\phi_{x}^{2} +
(\mathcal{H}_{x}\phi)^{2}\right)\right] + \frac{1}{2}
\phi_{t}^{2} \mathcal{H}_{xxx}\phi + \phi_{t} \phi_{tx}
\mathcal{H}_{xx}\phi \nonumber\\ &\quad\quad \ \ + \left. \phi_{x}
\partial_{x}\left[\phi_{t} \mathcal{H}_{tx}\phi -
\frac{1}{2} \left(\phi_{x}^{2} +
(\mathcal{H}_{x}\phi)^{2}\right)\right]\right\} +
{O}(\varepsilon^{3}) \quad \text{at}\ z=0,
\end{align}
which can be named as a Hilbert wave equation, general for small amplitude solutions and not restricted to narrow wavepackets. A version of this equation up to 
${O}(\varepsilon ^{1})$ appeared in Matsuno (Reference Matsuno1992), though deduced in a more involved fashion and without the goal to derive the NLS (1.11) for which 
${O}(\varepsilon ^{2})$-terms are required in (1.13). As known from perturbation theory, the higher-order nonlinear terms in (1.13) give rise to multiple scales in the problem in addition to the original 
$t$ and 
$x$ – this mechanism is different from the narrow wavepacket assumption causing 
$(T,\tau,X,Z)$-scales in the linear setting (§ 1.1). Equation (1.13) can be seen as a compact analogue of the Zakharov (Reference Zakharov1968) equation, but, as pointed by the author himself, triad resonances produce zero denominators and thus render his equation difficult to apply to capillary-gravity waves. Equation (1.13) and its straightforward generalization to the presence of surface tension, which would lead to a coupled system for 
$\phi (t,x)$ and 
$\eta (t,x)$, are not plagued with such difficulties. Finally, similar to the work of Dysthe (Reference Dysthe1979) on extending the perturbation analysis of the NLS (1.11) to the third order in 
$\varepsilon$, (1.13) readily admits the same extension.
Derivation of (1.11) from (1.13) becomes elementary since for the wavepacket solution
\begin{equation} \phi \sim e^{\mathrm{i} (x - t)} \left[\varPhi_{0}(T,\tau,X) + \varepsilon \varPhi_{1}(T,\tau,X) + \ldots\right] + \mathrm{c.c.} \equiv \phi_{0} + \varepsilon \phi_{1} + \cdots + \mathrm{c.c.}, \end{equation}
the HT acts on the narrow bandwidth waveform and not on slow modulation according to Bedrosian's theorem (King Reference King2009) because for 
${\rm e}^{\mathrm {i} x} \varPhi (X)$ the corresponding Fourier transforms 
$2 {\rm \pi}\delta (k-1)$ and 
$({1}/{\varepsilon } )\hat {\varPhi }({k}/{\varepsilon })$ of the involved factors are well-separated in the Fourier space, one being sharply peaked at 
$k=1$ as the Dirac's delta function 
$\delta(k-1)$ and the other having an 
${O}(\varepsilon )$-support around 
$k=0$. At the leading order, we find a ‘wave’ equation
\begin{equation} {O}(\varepsilon^{0}): \ \phi_{0 tt} + \mathcal{H}_{x}\phi_{0} = 0. \end{equation}
At the next order, keeping in mind that 
$\mathcal {H}$ acts on the fast scale 
$x$, we get
\begin{equation} {O}(\varepsilon^{1}): \ \phi_{1 tt} + \mathcal{H}_{x} \phi_{1} + \mathcal{H}\phi_{0 X} + 2 \phi_{0 tT} = \text{quadratic nonlinearity}, \end{equation}which entails solvability condition (1.5a). At the final order, where the NLS appears, we deduce
\begin{equation} {O}(\varepsilon^{2}): \ \phi_{2 tt} + \mathcal{H}_{x}\phi_{2} + \mathcal{H}\phi_{1 X} + 2 \phi_{1 tT} + \phi_{0 TT} + 2 \phi_{0 t \tau} = \text{nonlinearity}; \end{equation}
here we may put 
$\phi _{1}=0$ as per the discussion in § 1.1. As a result, the solvability condition on 
$\phi _{0}$ becomes 
$\varPhi _{0 TT} - 2 \mathrm {i} \varPhi _{0 \tau } = 0$, which with (1.5a) can be brought to the linear part of the NLS (1.5b). The nonlinearity in (1.11) then stems from the cubic terms in (1.13), which involve fast derivatives and the 
${O}(\varepsilon ^{0})$-solution only thereby making computations concise.
The linear wave equation (1.15) is obviously non-local, because free-surface perturbations propagate not only along the free surface, but also through the bulk. Its general solution is a combination of right- and left-travelling dispersive waves 
$\phi _{0}(t,x) = \int _{-\infty }^{+\infty }{\hat {\phi }_{0}(k) \exp [{\mathrm {i} (k x \mp \omega (k) t)}] \, \mathrm {d} k}$ and is illustrated in figure 1(a) comparatively with that of the ordinary wave equation 
$\phi _{0 tt} - \phi _{0 xx} = 0$: while the initial positive Cauchy pulse is breaking into two sign-definite hump-solutions translating without any change by the ordinary wave equation, the linear Hilbert wave equation (1.15) produces a faster spreading at long wavelengths with the waves travelling in the opposite directions constantly changing their shape because of dispersion – different wavenumbers travel with unequal speeds and hence the initial positivity of the sum of those at 
$t=0$ is not preserved. Numerical simulation of (1.13) with nonlinear terms and small 
$\varepsilon$ initially follows the solution of (1.15), but at some point exhibits a short-wave instability. The latter is easy to understand by ‘linearizing’ about the solution of the linear wave equation 
$\phi _{0}$, i.e. taking 
$\phi = \phi _{0} + \phi ^{\prime }$ we find for the growth rate 
$\lambda $ of the perturbation 
$\phi ^{\prime } \sim \exp ({\lambda t + \mathrm {i} k x})$ in the short-wave limit:
\begin{equation} \lambda^{2} \sim \varepsilon \phi_{0 t} k^{2}, \quad k \rightarrow \infty. \end{equation}
Therefore, if the free-surface acceleration 
$\phi _{0 t}$ is positive, one observes a Rayleigh–Taylor instability albeit modified by the free-surface curvature (Tanveer Reference Tanveer1991; Krechetnikov Reference Krechetnikov2009), which is short-wave and not regularized here due to the absence of surface tension, thus rendering it an ill-posed problem in the sense of Hadamard (Reference Hadamard1923). Hence, any application of (1.13) requires limiting the 
$k$-bandwidth from above (Joseph & Saut Reference Joseph and Saut1990) as in the presented derivation of the NLS (1.11).
The HT approach can also be developed for arbitrary depth h as envisaged by Matsuno (Reference Matsuno1992): imposing a no-penetration BC at 
$z = -h: \ \phi _{z}(t,x,-h) = 0$, the generalization of (1.12) is
\begin{equation} \phi_{z}(t,x,0) = \int_{-\infty}^{+\infty}{\phi_{x^{\prime}}(t,x^{\prime}) K^{h}(x-x^{\prime}) \, \mathrm{d} x^{\prime}} \equiv \mathcal{H}^{h} \phi_{x}(t,x) \end{equation}
with the kernel 
$K^{h}(x) = ({1}/{2 h}) \operatorname {csch}({{{\rm \pi} x}/{2 h}})$, which in the limit of infinite depth recovers the Cauchy kernel 
$K(x)$ in (1.12). To take the shallow-water limit, we integrate (1.19) by parts:
\begin{equation} z=0: \ \phi_{tt}(t,x) ={-} \phi_{z}(t,x,0) ={-} \int_{-\infty}^{+\infty}{\phi_{x^{\prime}x^{\prime}}(t,x^{\prime}) \delta^{h}(x-x^{\prime}) \, \mathrm{d} x^{\prime}}, \end{equation}
where 
$\delta ^{h}(x) = - ({1}/{{\rm \pi} h}) \ln |\tanh {({{\rm \pi} x}/{4 h})}|$ is a sequence approximating Dirac's delta function, though an unusual one being unbounded at 
$x=0$ even for finite 
$h$; however, 
$\int _{-\infty }^{+\infty }{\delta ^{h}(x) \, \mathrm {d} x} = 1$ for any 
$h$ and 
$\lim _{h \rightarrow 0}\int _{-\infty }^{+\infty }{f(x) \delta ^{h}(x) \, \mathrm {d} x} = f(0)$ for any test function 
$f(x)$ as required. Expanding in 
$h$, the linear part of (1.20) in the shallow-water limit becomes a local wave equation
\begin{equation} \phi_{tt} = \phi_{xx} + \tfrac{1}{3} h^{2} \phi_{xxxx} + \ldots, \end{equation}
i.e. the non-zero shallowness is responsible for weak dispersion. With the change of variables 
$(t,x) \rightarrow (\xi =x-t, \tau =h^{2} t)$, the above equation recovers the linear part of the Korteweg-de Vries (KdV) equation, 
$\phi _{\tau } + \frac {1}{6} \phi _{\xi \xi \xi } = \ldots$, and the nonlinear advection term arises from the quadratic nonlinearity in (1.13) with 
$\mathcal {H} \rightarrow \mathcal {H}^{h}$. Therefore, free surface itself is non-dispersive as per (1.21) and only non-zero bulk depth leading to non-locality is responsible for dispersion of water waves.
2. Axisymmetric water waves
In system (1.1) operators are now adapted to the cylindrical setting. Linearization of (1.2) results in the following form of the solution for the velocity potential analogous to (1.3):
\begin{align} \phi(t,r,z) &=
\int_{0}^{\infty}{\hat{\phi}_{0}(k) J_{0}(k r)
\exp({-\mathrm{i} \omega(k) t + k z}) \, k\, \mathrm{d} k}
+ \mathrm{c.c.} \nonumber\\ &\approx \varepsilon^{1/2}
\exp\left({\mathrm{i} (r - t) + z - \mathrm{i}
\frac{\rm \pi}{4}}\right) \frac{1}{\sqrt{2 {\rm \pi}R}}
\nonumber\\ &\quad \times
\int_{-\infty}^{\infty}\hat{\phi}_{0}(\kappa) k^{1/2}
\exp\left({\mathrm{i} \left[\kappa \left(R -
\omega_{0}^{\prime} T\right) -
\frac{\omega_{0}^{\prime\prime}}{2} \kappa^{2} \tau\right]
+ \kappa Z}\right) \nonumber\\ &\quad \times\left(1 -
\frac{\mathrm{i} \varepsilon}{8 k R} - \frac{9
\varepsilon^{2}}{128 k^{2} R^{2}} + \ldots\right)
\mathrm{d} \kappa + \mathrm{c.c.},
\end{align}
where on the right-hand side we showed only the outward-travelling wave as its contribution to the integral in (2.1) dominates at large distances as per the stationary-phase analysis (Koshlyakov, Smirnov & Gliner Reference Koshlyakov, Smirnov and Gliner1964), and 
$\hat {\phi }_{0}(k)$ is the Hankel transform of the IC 
$\phi _{0}(0,r,0)$; we also expanded the Bessel function 
$J_{0}(k r)$ of zero order away from the origin, 
$k r \gg 1$, used a narrow wavepacket 
$k = 1 + \varepsilon \kappa$ approximation, and introduced 
$R = \varepsilon r$. The solution (2.1) clearly illustrates the lack of sharp signals in two dimensions (2-D), i.e. does not satisfy the Huygens’ principle in the narrow sense (Hadamard Reference Hadamard1923), which in the context of waves on deep water is also enhanced by dispersion, the effect present in one dimension (1-D) as well (cf. figure 1a).
Equation (2.1) can be phrased as a formal expansion:
\begin{equation} \phi(t,r,z) =
\exp({\mathrm{i} (r - t) + z})
\left[\!\varPhi_{0} + \varepsilon
\varPhi_{1} + \varepsilon^{2}
\varPhi_{2} + \ldots\!\right](T,\tau,R,Z) + \mathrm{c.c.},
\end{equation}
i.e. 
$\phi (t,r,z)$ is a travelling wave 
$\exp [{\mathrm {i} (r - t)}]$ modulated with envelope amplitudes 
$\varPhi _{i}(T,\tau,R,Z)$, 
$i=0,1,\ldots$, each travelling with the same group speed 
$\omega _{0}^{\prime }$ as per the exponential in the integrand of (2.1). Identifying the first envelope function 
$\varPhi _{0}(T,\tau,R)=\varPhi _{0}(T,\tau,R,0)$ from (2.1) we find by direct differentiation that it satisfies the condition of a wavepacket propagating at the group velocity and the linear part of the envelope equation in analogy to (1.5) at 
$Z=0$:
\begin{equation} \varPhi_{0T} + \omega_{0}^{\prime} \left(\varPhi_{0R} + \frac{1}{2 R} \varPhi_{0}\right) = 0, \quad \varPhi_{0\tau} - \frac{\mathrm{i} \omega_{0}^{\prime\prime}}{2} \left(\varPhi_{0RR} + \frac{1}{R} \varPhi_{0R} - \frac{1}{4 R^{2}} \varPhi_{0}\right) = 0, \end{equation}
respectively. With the transformation 
$\varPhi _{0}(T,\tau,R) = R^{-1/2} \tilde {\varPhi }_{0}(T,\tau,R)$ the above equations reduce to the familiar translationally invariant plane counterparts (1.5). Thus, the potential 
$- {1}/{4 R^{2}}$ is required to guarantee mass conservation 
$\int _{0}^{\infty }{\eta _{0}(r) r \,\mathrm {d} r} = \int _{0}^{\infty }{\eta (t,r) r \,\mathrm {d} r}$ via modifying the cylindrical wave amplitude with the weight 
$R^{-1/2}$ as it travels to or from the origin; we call such a translational invariance weighted. The presence of the potential in (2.3b) implies that the NLS is not covariant in view of being a reduced-order model rather than a fundamental physical law (Krechetnikov Reference Krechetnikov2024).
2.1. Indeterminacy of the multiple-scales approach
As we recognized from (2.1), a weakly nonlinear model, if exists, must reside at the radial distances 
$R = \varepsilon r = {O}(1)$. As a result, at the leading order 
${O}(\varepsilon ^{0})$ we get the system defining the carrier wave identical to that in the plane case (1.6) with 
$x \rightarrow r$ and the outward-travelling wavepacket solution being 
$\phi _{0} = \varPhi _{0}(T,\tau,R,Z) \exp ({\mathrm {i} (r - t) + z}) + \mathrm {c.c.}$. At the next order, some differences occur:
\begin{equation} {O}(\varepsilon^{1}):
\begin{cases} \phi_{1 zz} + \phi_{1 rr} ={-} 2
\left(\phi_{0 zZ} + \phi_{0 rR}\right) - \dfrac{\phi_{0
r}}{R} & \text{for} \ z \leqslant 0, \\ \phi_{1 z} +
\phi_{1 tt} = \,\,\, {-} 2 \phi_{0 t T} - \phi_{0 Z} - \eta_{0 t}
\phi_{0 z t} - \phi_{0 z} \phi_{0 z t} & \\ \quad\qquad\qquad\ \ \,\, -
\eta_{0} \left(\phi_{0 z tt} + \phi_{0 zz}\right) + \eta_{0
r} \phi_{0 r} - \phi_{0 r} \phi_{0 r t} & \text{at} \ z =
0. \end{cases} \end{equation}
The right-hand side of the Poisson equation in (2.4) brings about secular terms with 
$\exp [{\pm \mathrm {i} (r - t)}]$, the factors of which vanish provided that the no-resonance condition on 
$\varPhi _{0}(T,\tau,R,Z)$, analogous to (1.8), holds:
\begin{equation} \varPhi_{0 Z} + \mathrm{i} \left(\varPhi_{0 R} + \frac{\varPhi_{0}}{2R}\right) = 0 \quad \text{for} \ Z \leqslant 0; \end{equation}
differential consequence 
$\varPhi _{0 ZZ} + ({1}/{R})\partial _{R} (R \varPhi _{0 R}) = ({1}/{4 R^{2}})\varPhi _{0}$ of (2.5) shows that 
$\varPhi _{0}$ is no longer a harmonic function in the 
$(R,Z)$-space in contrast to its plane equivalent, which results in the loss of translational invariance and emergence of the 
$-{1}/{4 R^{2}}$-potential in (2.3b).
The BCs in (2.4) give rise to the no-resonance requirement (2.3a), which with (2.5) at 
$Z=0$ yields
\begin{equation} Z=0: \ \mathrm{i} \omega_{0}^{\prime} \varPhi_{0 Z} + \varPhi_{0 T} = 0. \end{equation}The solvability condition (2.3a) can be easily integrated to produce
\begin{equation} \varPhi_{0}(T,\tau,R) = \frac{1}{R^{1/2}} \tilde{\varPhi}_{0}\left(\tau,\xi\right), \ \text{with} \ \xi = R-\omega_{0}^{\prime} T, \end{equation}
so that 
$\tilde {\varPhi }_{0}$ satisfies 
$\tilde {\varPhi }_{0 T} + \omega _{0}^{\prime } \tilde {\varPhi }_{0 R} = 0$, i.e. the equivalent of (1.5a), and shows that the first two terms in (2.3a) represent advection accounting for the envelope travelling at the group velocity and the last term is responsible for dilution decreasing the wavepacket amplitude with the radial distance as 
$\varPhi _{0} \sim R^{-1/2}$. Similar to the plane case (§ 1.1), the solution 
$\phi _{1} = \varPhi _{1}(T,\tau,R,Z) \exp ({\mathrm {i} (r - t) + z}) + \mathrm {c.c.}$ is the same as for 
$\phi _{0}$ with no other harmonics, 
$\exp [{\mathrm {i} n(r-t)}]$, 
$n=0,\pm 2,\ldots$, excited.
Finally, the Poisson equation at the order at which an envelope equation appears,
\begin{equation} {O}(\varepsilon^{2}): \ \phi_{2 zz} + \phi_{2 rr} ={-} 2 \left(\phi_{1 zZ} + \phi_{1 rR}\right) - \left(\phi_{0 ZZ} + \phi_{0 RR}\right) - \frac{\phi_{1 r}}{R} - \frac{\phi_{0 R}}{R} \quad \text{for} \ z \leqslant 0, \end{equation}
leads to the solvability condition on 
$\varPhi _{1}$, simplified with the differential consequence of (2.5),
\begin{equation} \varPhi_{1 Z} + \mathrm{i} \left(\varPhi_{1 R} + \frac{1}{2 R} \varPhi_{1}\right) ={-} \frac{1}{8 R^{2}} \varPhi_{0} \quad \text{for} \ Z \leqslant 0. \end{equation}
Hence, despite that the form of the solution for 
$\phi _{1}$ is the same as for 
$\phi _{0}$, which in the plane case allowed us to put 
$\phi _{1}=0$, (2.9) prevents one from doing that in view of entanglement of 
$\phi _{0}$ and 
$\phi _{1}$ owing to the cylindrical geometry. A posteriori, we realize that the term 
${O}(R^{-1})$ in (2.1) must contribute to the equation for 
$\varPhi _{1}$ by coupling to 
$\varPhi _{0}$. The combined BC at this order is
\begin{align} {O}(\varepsilon^{2}): \
\phi_{2 z} &+ \phi_{2 tt} ={-} 2 \phi_{1 t T} - 2 \phi_{0 t
\tau} - \phi_{1 Z} - \phi_{1 TT}\nonumber\\ &= {\rm
e}^{\mathrm{i} (r - t)} \left(2 \mathrm{i} \varPhi_{1 T} -
\varPhi_{1 Z}\right) + {\rm e}^{\mathrm{i} (r - t)} \left(2
\mathrm{i} \varPhi_{0 \tau} - \varPhi_{0 TT}\right) +
\mathrm{c.c.} \quad \text{at} \ z=0.
\end{align}
The crux of the matter is that the solvability condition ensuing from the right-hand side of (2.10) admits two interpretations, i.e. either to treat 
$\varPhi _{0}$ and 
$\varPhi _{1}$ in a coupled way or independently, which is indeterminacy of the multiple-scale approach revealed here in the axisymmetric geometry. Even though (2.9) suggests otherwise, if one chooses to comply with (2.1), in which the wavepacket width 
$\varepsilon$ is inversely proportional to the radial distance scale, 
$r \sim \varepsilon ^{-1}$, it seems that only one choice is valid as all the relevant information is already encapsulated in that linear solution (2.1). The latter imposes the condition 
$\mathrm {i} \omega _{0}^{\prime } \varPhi _{1 Z} + \varPhi _{1 T} = 0$, same as (2.6) on 
$\varPhi _{0}$, thus nullifying the first brackets on the right-hand side of (2.10) and leading to
\begin{equation} \varPhi_{1 T} + \omega_{0}^{\prime} \left(\varPhi_{1 R} + \frac{\varPhi_{1}}{2 R}\right) = \frac{\mathrm{i} \omega_{0}^{\prime}}{8 R^{2}} \varPhi_{0} \quad \text{at} \ Z=0, \end{equation}
whereas the second brackets in (2.10) produce (2.3b). While the left-hand side of (2.11) admits travelling solutions as at the previous order (2.3a), consistent with the fact that there exists only one physical mechanism (and hence speed) for wave propagation, (2.11) also has a particular solution. Indeed, with 
$\varPhi _{1} = R^{-1/2}\tilde {\varPhi }_{1}(T,\tau,R)$ and (2.7) from (2.11) we get
\begin{equation} \tilde{\varPhi}_{1 T} + \omega_{0}^{\prime} \tilde{\varPhi}_{1 R} = \frac{\mathrm{i} \omega_{0}^{\prime}}{8 R^{2}} \tilde{\varPhi}_{0} \quad \text{at} \ Z=0, \end{equation}
which admits a particular solution 
$\tilde {\varPhi }_{1} = - ({\mathrm {i}}/{8 R}) \tilde {\varPhi }_{0}(\xi )$, also propagating at the group velocity. Compared with the 
$\varPhi _{0}$-solution, the phase of 
$\varPhi _{1}$ alters by 
$-{{\rm \pi} }/{2}$ and amplitude decreases faster with 
$R$, which is consistent with (2.1) under closer consideration. Therefore, from this interpretation we infer the crucial difference from the plane case: the next-order amplitude 
$\varPhi _{1}$ is generated at the linear level by the coupling term 
$({1}/{R})\phi _{1 r}$ in (2.8), while in the plane case it would be generated by nonlinear interactions only leading to the clustered mode-structure of the solution.
The complete nonlinear version of the combined boundary condition (2.10) (not shown due to excessive number of terms), simplified with (2.5), (2.3a), (2.9) and their differential consequences, gives rise to the no-resonance condition in the form of the Gross–Pitaevskii equation:
\begin{equation} - 2 \mathrm{i} \varPhi_{0 \tau} + \frac{1}{4} \left[\varPhi_{0 RR} + \frac{\varPhi_{0 R}}{R} - \frac{1}{4} \frac{\varPhi_{0}}{R^{2}}\right] + 4 |\varPhi_{0}|^{2} \varPhi_{0} = 0 \quad \text{at} \ Z = 0, \end{equation}
which, compared with the 2-D NLS originally derived in Cartesian coordinates and then rewritten in cylindrical ones, is amended with an inverse-square potential 
$- {1}/{4 R^{2}}$. The weighted translational invariance no longer takes place because of the nonlinearity, i.e. naive substitution of (2.7) into (2.13) brings up incompatibility with the dependence 
$\tilde {\varPhi }_{0}(\xi )$ and thus the region of applicability of (2.13) in the 
$R$-space shrinks to zero if one seeks travelling modulated waves. Thus, it is the order 
${O}(\varepsilon ^{2})$ where the weighted translational invariance (2.7) is lost, which is expected as the balance of nonlinearity and dispersion can happen only for 
$R={O}(1)$ akin to the derivation of the concentric KdV (cKdV) equation (Iordansky Reference Iordansky1959; Maxon & Viecelli Reference Maxon and Viecelli1974), i.e. the resulting weakly nonlinear model cannot be uniformly valid everywhere (for all 
$R$) compared with that in the Cartesian setting.
However, we recall that multiple scales are, in fact, related so that the 
$\xi$ variable introduced in (2.7) can be rewritten as 
$R = \xi + ({1}/{2 \varepsilon } )\tau$ similar to what is done in the derivation of the cKdV (Johnson Reference Johnson1980). Hence, substituting (2.7) into (2.13) as well as rescaling 
$\tilde {\varPhi }_{0} \rightarrow \varepsilon ^{-1/2} \tilde {\varPhi }_{0}$ yields the concentric NLS (cNLS) with a 
$\tau$-dependent factor:
\begin{equation} \mathrm{i} \tilde{\varPhi}_{0 \tau} - \frac{1}{8} \tilde{\varPhi}_{0 \xi\xi} - \frac{2}{\tau} |\tilde{\varPhi}_{0}|^{2} \tilde{\varPhi}_{0} = 0 \quad \text{at} \ Z = 0, \end{equation}
where the effect of the 
$-{1}/{4 R^{2}}$-potential disappears and 
$\xi \in (-\infty,+\infty )$. Because the wavepacket travels at the group velocity, thus necessarily operating on the slow scale 
$R$ as opposed to cKdV, in which the 
$\xi$-frame is based on the fast scale 
$r$, there is no R-version of the 
$\tau$-equation (2.14) compared with the shallow-water equivalent (Johnson Reference Johnson1980). Same as the NLS (1.11), under the appropriate decay conditions for 
$\xi \rightarrow \pm \infty$ the cNLS possesses the conservation law 
$\mathrm {d} N / \mathrm {d} \tau \equiv ({\mathrm {d}}/{\mathrm {d} \tau })\int {|\tilde {\varPhi }_{0}|^{2} \, \mathrm {d} \xi } = 0$ for the number of particles 
$N$ named so in analogy to quantum mechanics, but the energy conservation reads 
$\mathrm {d} H / \mathrm {d} \tau = - ({1}/{\tau ^{2}})\int {|\tilde {\varPhi }_{0}|^{4} \, \mathrm {d} \xi }$ with the time-dependent ‘Hamiltonian’ 
$H = - \frac {1}{8}\int _{-\infty }^{+\infty }{|\tilde {\varPhi }_{0}|^{2} \, \mathrm {d} \xi } + ({1}/{\tau } )\int _{-\infty }^{+\infty }{|\tilde {\varPhi }_{0}|^{4} \, \mathrm {d} \xi }$ reflecting the time-transient nature of the solution 
$\tilde {\varPhi }_{0}$.
Whereas (1.11) admits standing-wave (Stokes) solutions 
$\varPhi _{0}(\tau,\xi ) = \mu ^{1/2} {\rm e}^{-\mathrm {i} \mu \tau } \operatorname {sech}{[(8 \mu )^{1/2} \xi ]}$ which are parameterized with the propagation constant 
$\mu$ and generating travelling-wave solutions 
$\varPhi _{0}(\tau,\xi ) \rightarrow \exp [{- 4 \mathrm {i} V (\tilde {\xi } + {\tau V}/{2})}] \varPhi _{0}(\tau,\tilde {\xi })$ owing to translational symmetry 
$(\tau,\xi ) \rightarrow (\tau,\tilde {\xi } = \xi - V \tau )$, the cNLS (2.14) does not. Therefore, we can compare them only on a set of self-similar solutions. The NLS admits self-similarity 
$\varPhi _{0} = \tau ^{-1/2} \varPi (\zeta )$ with 
$\zeta = \tau ^{-1/2} \xi$, and cNLS – 
$\tilde {\varPhi }_{0} = \varPi (\zeta )$ with 
$\zeta = \tau ^{-1/2} \xi$, thus leading to
\begin{equation} - \frac{\mathrm{i}}{2} \left(\underline{\varPi} + \zeta \varPi^{\prime}\right) - \frac{1}{8} \varPi^{\prime\prime} - 2 |\varPi|^{2} \varPi = 0, \end{equation}
where the underlined term is present only in the NLS case; the primes stand for the derivatives w.r.t. 
$\zeta$. Expressing 
$\varPi = A {\rm e}^{\mathrm {i} \varphi }$, we get a system for the amplitude 
$A(\zeta )$ and phase 
$\varphi (\zeta )$:
$$\begin{align}\mathrm{Re}&:\quad A^{\prime\prime} - A \varphi^{\prime 2} - 4 \zeta A \varphi^{\prime} + 16 A^{3} = 0, \end{align}$$
$$\begin{align}\mathrm{Im}&:\quad A \varphi^{\prime\prime} + 2 A^{\prime} \varphi^{\prime} + 4 \zeta A^{\prime} + \underline{4 A} = 0. \end{align}$$
Both cases allow for odd and even solutions 
$A(\zeta )$ w.r.t. 
$\zeta \rightarrow -\zeta$ with the asymptotics
\begin{equation} \zeta \rightarrow +\infty : A \sim \left(\frac{\epsilon}{\zeta}\right)^{1/2} \cos{\varOmega} - \frac{3}{8}\!\left(\frac{\epsilon}{\zeta}\right)^{3/2} \sin{\varOmega} + \ldots, \ \varphi ={-} \zeta^{2} \mp \frac{1}{4 \zeta^{2}} \ln{|\sin{\varOmega}|} + \ldots, \end{equation}
and similar for 
$\zeta \rightarrow -\infty$; here 
$\epsilon = \sqrt {2}-1$, 
$\varOmega = \zeta ^{2} + {{\rm \pi} }/{4}$ and the only difference between the NLS and cNLS appears in the phase 
$\phi$ with the upper sign choice corresponding to the cNLS and the lower to the NLS. Equation (2.16) also possesses dispersive shock-wave (SW) solutions (Gurevich & Pitaevskii Reference Gurevich and Pitaevskii1974) potentially relevant to undular bores, for example, as shown in figure 1(b) – they are virtually indistinguishable between the NLS and cNLS cases. The cNLS (2.14) and the above analysis of self-similar solutions to (1.11), (2.14) provide a deep-water counterpart to an analogous study of self-similar solutions of KdV and cKdV (Johnson Reference Johnson1980) for shallow water, which also exhibit asymptotics that are slowly decaying and oscillatory at infinity.
2.2. Analogue of the HT in the axisymmetric case
As in the plane case (1.12), the idea is to express 
$\phi _{z}(t,r,0)$ in terms of 
$\phi (t,r)$ or 
$\phi _{r}(t,r)$ with the former choice proving to be more convenient in the axisymmetric case. To this end, from the solution 
$\phi (t,r,z) = \int _{0}^{+\infty }{k J_{0}(k r) \hat {\phi }(t,k) {\rm e}^{k z} \, \mathrm {d} k}$ of the Laplace equation (1.1a) in the lower half-space supplied with the Dirichlet BC 
$\phi (t,r,0) = \phi (t,r)$, we may calculate
\begin{equation} \phi_{z}(t,r) = \int_{0}^{+\infty}{k^{2} J_{0}(k r) \hat{\phi}(t,k) \, \mathrm{d} k}. \end{equation}
To simplify (2.18), let us resort to the convolution in the 
$\boldsymbol {r}=(r,\theta )$-plane (Baddour Reference Baddour2009):
\begin{align} \gamma(\boldsymbol{r}) &= \alpha(\boldsymbol{r}) ** \, \beta(\boldsymbol{r}) \equiv \iint_{\Bbb{R}^{2}}{\beta(\boldsymbol{r}-\boldsymbol{r}^{\prime}) \alpha(\boldsymbol{r}^{\prime}) \, \mathrm{d} \boldsymbol{r}^{\prime}} \nonumber\\ &= \int_{0}^{\infty}{\alpha(r^{\prime}) \, G(r-r^{\prime}) r^{\prime} \, \mathrm{d} r^{\prime}} = 2 {\rm \pi}\int_{0}^{+\infty}{\hat{\beta}(k) \hat{\alpha}(k) J_{0}(k r) k \, \mathrm{d} k}, \end{align}
where 
$G(r-r^{\prime }) = \int _{0}^{2{\rm \pi} }{\beta (\boldsymbol {r}-\boldsymbol {r}^{\prime }) \, \mathrm {d} \theta ^{\prime }}$, and 
$\hat {\beta }(k)$, 
$\hat {\alpha }(k)$ are the respective Hankel transforms of 
$\beta (r)$,
$\alpha (r)$, so that in the transformed space 
$\hat {\gamma }(k) = 2 {\rm \pi}\hat {\alpha }(k) \hat {\beta }(k)$. Applying this to (2.18) produces
\begin{equation} \phi_{z}(t,r) = \int_{0}^{+\infty}{\phi(t,r^{\prime}) G(r-r^{\prime}) r^{\prime} \mathrm{d} r^{\prime}}, \end{equation}
where the inverse transform of 
$\hat {\beta }(k)$ is 
$\beta (r) = {1}/{2 {\rm \pi}r^{3}}$ and, correspondingly,
\begin{equation} G(r-r^{\prime}) =
\int_{0}^{2{\rm \pi}}{\beta(\boldsymbol{r}-\boldsymbol{r}^{\prime})
\, \mathrm{d} \theta^{\prime}} ={-} \frac{1}{\rm \pi}
\frac{\rho(r,r^{\prime})}{(r-r^{\prime})^{2}
r^{\prime}}, \quad \rho(r,r^{\prime}) = \frac{2
r^{\prime}}{r+r^{\prime}} E\left(\frac{2 \sqrt{r
r^{\prime}}}{r +
r^{\prime}}\right);
\end{equation}
here 
$E({\cdot })$ is the complete elliptic integral of the second kind.
Now, we are in a position to recover the envelope equation (2.13). For brevity, we will focus on the linear terms since they are the source of indeterminacy in the multiple-scale approach (§ 2.1). The linear wave equation resulting from the BCs is
\begin{equation} z=0: \ \phi_{tt} ={-} \phi_{z}(t,r) ={-} \int_{0}^{+\infty}{\phi(t,r^{\prime}) G(r-r^{\prime}) r^{\prime} \, \mathrm{d} r^{\prime}}. \end{equation}
It transpires that the kernel 
$G(r-r^{\prime })$ should be almost 1-D for 
$r = {O}(\varepsilon ^{-1})$. In general, for any 
$r$, 
$r^{\prime }$ one may separate the singular and regular parts of the kernel (cf. figure 2a):
\begin{equation} G(r-r^{\prime}) ={-} \frac{1}{{\rm \pi} r^{\prime}} \frac{1}{(r-r^{\prime})^{2}} + \frac{1}{2 {\rm \pi}r r^{\prime}} \frac{1}{r-r^{\prime}} + \frac{1}{8 {\rm \pi}r^{2} r^{\prime}} \ln{|r-r^{\prime}|} + \text{regular part}, \end{equation}where the singular part yields
\begin{align}
\int_{0}^{+\infty}\phi(t,r^{\prime}) G(r&-r^{\prime})
r^{\prime} \, \mathrm{d} r^{\prime} ={-} \frac{1}{\rm \pi}
\int_{0}^{+\infty}{\frac{\phi(t,r^{\prime})}{(r-r^{\prime})^{2}}
\, \mathrm{d} r^{\prime}} + \frac{1}{2 {\rm \pi}r}
\int_{0}^{+\infty}{\frac{\phi(t,r^{\prime})}{r-r^{\prime}}
\, \mathrm{d} r^{\prime}} + \mathrm{h.o.t.} \nonumber\\ &=
\frac{\partial}{\partial r} \mathcal{H}\phi(t,r) +
\frac{\varepsilon}{2 R} \mathcal{H}\phi(t,r) +
\frac{\varepsilon^{2}}{8 R^{2}} \partial_{r}^{-1}
\mathcal{H}\phi(t,r) + \mathrm{h.o.t.};
\end{align}
here the lower limit of integration can be changed to 
$-\infty$ because the modulated solution is non-zero only far from the origin, 
$r = {O}(\varepsilon ^{-1})$, and inessential for our discussion terms are of higher order (h.o.t.). Next, applying the introduced earlier multiple-scale ansatz and expanding the solution in series, at the first three orders we obtain, in analogy to (1.15)–(1.17) in the plane case, the solvability conditions
$$\begin{align} {O}(\varepsilon^{0})&: \ \phi_{0 tt} + \mathcal{H}_{r}\phi_{0} = 0, \end{align}$$
$$\begin{align}{O}(\varepsilon^{1})&: \ 2 \phi_{0 t T} + \left(\mathcal{H}\phi_{0 R} + \frac{1}{2 R}\mathcal{H}\phi_{0}\right) = 0, \end{align}$$
$$\begin{align}{O}(\varepsilon^{2})&: \ 2 \phi_{0 t\tau} + \phi_{0 TT} = 0, \end{align}$$
i.e. the same as in § 2.1 after we removed the indeterminacy from the 
${O}(\varepsilon^{2})$ condition with the help of (2.11) emerged here due to the second and third terms in (2.23); in (2.25c) one may replace 
$\phi _{0 TT}$ with the differential consequence of the 
${O}(\varepsilon ^{1})$-equation and represent 
$\phi _{0} = \varPhi _{0}(T,\tau,R) {\rm e}^{\mathrm {i} (r-t)} + \mathrm {c.c.}$ thus furnishing the linear part of (2.13). Therefore, it is the second and third terms in (2.23) that are absent in the plane case and responsible for coupling of 
$\phi _{0}$, 
$\phi _{1}$ and hence for the indeterminacy in the axisymmetric case we encountered in § 2.1.
(a) Regular part of the kernel (2.23) (solid) vs the singular behaviour of the kernel itself (dashed) – here 
$G(r/r^{\prime }) = r^{\prime 3} G(r-r^{\prime })$; (b) the HT of the derivative of the Cauchy pulse (dashed) in the plane case (blue, 
$c=x$) vs axisymmetric case (red, 
$c=r$).
To elucidate this difference from the plane analogue (1.12), we integrate (2.20) by parts:
\begin{equation} \phi_{z}(t,r) = \frac{1}{\rm \pi} \int_{0}^{+\infty}{\frac{1}{r-r^{\prime}} \frac{\mathrm{d}}{\mathrm{d} r^{\prime}}\left[\rho(r,r^{\prime}) \phi(t,r^{\prime})\right] \mathrm{d} r^{\prime}}, \end{equation}
where we took into account that 
$\rho(r,0)=0$. Equation (2.26) is the counterpart of (1.12) and its form explains why starting with the Dirichlet data 
$\phi (t,r)$ rather than 
$\phi _{r}(t,r)$ is a more convenient choice in the axisymmetric case. Hence, if we are to apply the plane HT to the derivative of the Cauchy pulse 
$\phi (x) = {1}/({1+x^{2}})$ for comparison, we should act with (2.26) on (the half of) the Cauchy pulse 
$\phi (r) = {1}/({1+r^{2}})$, 
$r>0$. Compared with the kernel in (1.12) 
$\sim 1/{\rm \pi}(x-x')^2$ in the plane case, the kernel (2.21a,b) in addition to the same behaviour 
$\sim 1/{\rm \pi} r' (r-r')^2$ in the neighbourhood of the singularity 
$r=r^{\prime }$, where the factor 
$1/r^{\prime }$ cancels out with the cylindrical measure in (2.20), contains extra singular terms per (2.23) as well as a non-zero regular part. Hence, convolution is only weakly localized in the cylindrical case, provoking coupling between the leading 
$\phi _{0}$ and next-order 
$\phi _{1}$ effects we observed in § 2.1 and, as a result, radial spreading of the wave relative to the plane case (cf. figure 2b).
3. Conclusions
In water-wave mechanics the HT has been previously encountered in internal waves (Davis & Acrivos Reference Davis and Acrivos1967) and traditionally used to analyse experimental water-wave data (Melville Reference Melville1983; Huang, Shen & Long Reference Huang, Shen and Long1999). As we demonstrated both in the plane and axisymmetric cases, the advantage of the HT approach is also in significant reduction of algebra in the course of derivation of envelope equations for surface waves (§ 1.1and § 2.1), in particular, in the number of solvability conditions, because one does not need to handle the solutions in the bulk as in the now classical multiple-scale approach. Moreover, in the axisymmetric case (§ 2.2) it highlights the source of the indeterminacy peculiar to this geometric setting (§ 2.1). Also, the HT technique produces a general bandwidth equation (1.13), which has certain advantages compared with Zakharov's equation, in particular, avoiding resonances and multi-dimensional Fourier integrals, and therefore is worth further exploration, for example, in the context of wave turbulence (Nazarenko & Lukaschuk Reference Nazarenko and Lukaschuk2016). Generalization of the HT to a finite depth (§ 1.2) elucidates how one transitions from a purely dispersive and non-local wave equation (1.13) on deep water to a weakly dispersive and local wave equation (1.21) on shallow water. Finally, the HT analogy in the axisymmetric case (§ 2.2) is instrumental in getting an insight into coupling of 
$\phi _{0}$- and 
$\phi _{1}$-potentials, the effect absent in the plane case, and leads to the cNLS (2.14). Extension of the presented approach to general 2-D waves can be based on a generalization of (1.12) to the 
$(x,y)$-plane:
\begin{equation} \phi_{z}(x,y,0) ={-} \frac{1}{2 {\rm \pi}} \iint_{-\infty}^{+\infty}{\phi(x^{\prime},y^{\prime},0) \left[(x-x^{\prime})^{2} + (y-y^{\prime})^{2}\right]^{{-}3/2} \,\mathrm{d} x^{\prime} \,\mathrm{d} y^{\prime}}, \end{equation}which in the radial symmetry reduces to (2.20), but otherwise requires one to analyse a Poisson-type integral transform leading to the 2-D Schrödinger equation (Zakharov Reference Zakharov1968).
Acknowledgements
In memory of Professor G.M. Bud Homsy, who kindly read an earlier version of the manuscript. The author is also grateful to Professor T. Colonius, Dr A. Zelnikov and Q. Farr for the discussions on the subject of this study.
Funding
The work was supported by the Natural Sciences and Engineering Research Council of Canada under grant no. 04374.
Declaration of interests
The author reports no conflict of interest.
