The classic problem of ship waves of infinitely small amplitudes is studied within the recently developed approach of reference solutions [V. G. Gnevyshev and S. I. Badulin, Moscow Univ. Phys. Bull. 72, 415 (2017)]. The evolution of narrow-banded Gaussian wavetrains of a finite volume is considered as an alternative to the conventional inherently point-wise tracking of wavetrains within the asymptotic methods of the stationary phase or the steepest descent. The approach allows for avoiding general problems of the methods: occurrence of singularities of wave fields. The non-singular solutions for a stationary ship wake are presented in an analytical form in terms of two key dimensionless quantities: the Froude number and the aspect ratio of the initial (boundary) domain of wave generation. Systems of transverse and diverging ship waves as well as amplitude and phase effects of their transformation can be analyzed separately within this approach for small Froude numbers. The essential role of dispersion of three-dimensional water waves is emphasized and detailed.

This paper presents asymptotic solutions to the classic problem of stationary deep water ship waves first formulated by Lord Kelvin.1 The Kelvin analysis of wave kinematics predicts universal wake patterns in the form of arms of chevron accompanied by a system of transverse arcs inside a gusset of total opening at approximately 39° (twice the Kelvin angle, θK ≈ 19.47°). This analysis was later extended2–5 by studies of amplitude distributions based on asymptotical approaches of the stationary phase6 (see also the work of Thomson7) and the steepest descent.8,9 It should be mentioned that these approaches showed their efficiency in many wave problems, mostly one-dimensional. An extension of these methods to two and more dimensions is painstaking or even impossible. Quite often, the solution of the multi-dimensional problem is replaced by a unidirectional one to avoid serious mathematical difficulties.

Wave caustics are the main difficulty of the stationary phase method (as well as of its counterpart, the steepest descent approach). The difficulty lies in an essentially point-wise nature of the method when the evolution of a spatial wave pattern of finite volume is replaced by one of the characteristic points that provides phase stationarity. In this way, the concept of wave group velocity appears6,10 when the evolution can be adequately described as a propagation of the characteristic point of the wavetrain at the group velocity. The evolution of wave amplitude is accounted for by the dispersion of group velocity (second derivatives of the dispersion relationship). The spatial distribution of the wavetrain remains beyond this method. When wave dispersion vanishes at a caustic (or at its projection when reduced to an omnidirectional problem), this approach faces a mathematical problem: there appears singularity of wave amplitude which can be removed by accounting for higher-order effects of wave dispersion. A physical way of avoiding the “mathematical infinity” “to blur or to smooth it down” was realized by Lord Kelvin1 (p. 425) quite well by using explicit wavetrain distributions2 [Eq. (95)]. This physical idea found its mathematical treatment in many works on asymptotic methods and perturbation theory (Ref. 11 and references therein).

The novel approach recently sketched by the authors,12 the so-called reference solution method, follows the idea of Lord Kelvin. According to this approach, the reference wavetrains are Gaussian, which allow one to develop an analytic approach to many problems of wave propagation in multidimensional dispersive media. The resulting asymptotic solutions are shown to be finite and continuous in the whole space.

In Sec. II, we outline features of the method of reference solutions for dispersive waves of infinitely small amplitude. The general formulas will be used in Sec. III for the classic problem of stationary ship waves. Solutions will be presented in an explicit form where particular terms have transparent physical meaning. In this way, the effects of amplitude and phase transformation can be analyzed separately. Graphical illustrations unveil the general effects of ship wake formation. In conclusion, we discuss the prospects of the proposed approach.

The new approach to the analysis of linear wave packets in the case of two or more spatial dimensions has been proposed recently by the authors.12 In this section, we present the basics of the method and general formulas for wavetrains in dispersive multidimensional media.

Following the study of Gnevyshev and Badulin,12 consider solutions for waves of infinitely small amplitude in the form of linear superposition of harmonics in two-dimensional media (e.g., water waves),

ηx,t=Re(2π)1F(k)exp[ikxΩ(k)t]dk.
(1)

Here, η(x, t) is a wave variable (e.g., surface elevation of water waves), Re is the real part of the corresponding expression, and F(k) is the Fourier image of the function η(x, t) at initial time t = 0. The wavevector k = (k, l) and frequency Ω are assumed to obey the dispersion relation Ω(k).

Consider the initial distribution f(x) = η(x, t = 0) in (1) in the form of a Gaussian wavetrain,

f(x,y)=Re1(ΔxΔy)expik0x+il0yx22(Δx)2y22(Δy)2.
(2)

Note that at the limit Δx → 0, Δy → 0, (2) is reduced to a δ-like function,

limΔx0,Δy0f(x,y)=δ(xx0)δ(yy0)2cos(k0x+l0y).
(3)

This property will be used essentially to construct solutions for the ship wake of a finite bandwidth. In the k-space distribution, (2) is a normalized Gaussian pulse,

F(k,l)=expκx22(Δk)2κy22(Δl)2,
(4)

where the wavevector κ = (κx, κy) denotes the deviation from the carrier wavevector k0 = (k0, l0), i.e., κx = (kk0), κy = (ll0), and the spectral bandwidth is related to its dimension in a coordinate space as follows:

Δk=(Δx)1;Δl=(Δy)1.
(5)

Assuming the pulse in (4) is narrow-banded, i.e., Δk/k0 ≪ 1, Δl/l0 ≪ 1, one can approximate the dispersion relation Ω(k) in (1) by its Taylor series in the vicinity of the wave packet carrier harmonic k0 (Ω0=Ωk=k0),

Ω(k)=Ω0+Ωkk=k0κx+Ωlk=k0κy+12Ωkkk=k0κx2+Ωllk=k0κy2+2Ωklk=k0κxκy+O(|κ|3).
(6)

Herein, we use the following notations for the second order differentials associated with the bandwidth of the reference function (4):

μx=Ωkkk=k0(Δk)2;μy=Ωllk=k0(Δl)2;μxy=Ωklk=k0(ΔkΔl).
(7)

The resulting approximation for (1) takes a form of a modulated wavetrain,

η(x,t)=Re+expBλ,λ2iλ,ξdλ×expi𝜗0,
(8)

with the phase of the carrier wave being

𝜗0=k0xΩ0t,
(9)

and the dimensionless wavevector λ and the coordinate ξ that characterizes formally small deviations from the carrier harmonic k0 and its trajectory in a coordinate space being

λ=(κx/Δk,κy/Δl),
(10)

and

ξ=Δk(xΩkk=k0t),ΔlyΩlk=k0t.
(11)

The modulation function given by the Fourier integral in (8) is determined by the quadratic form with a dimensionless matrix,

B=1+itμxitμxyitμxy1+itμy.
(12)

The quadratic form exponent in (8) can be integrated explicitly8 as follows:

+expBλ,λ2iλ,ξdλ=2πdet(B)1/2expB1ξ,ξ2.
(13)

The real part of the matrix B is strictly positively determined that ensures convergence of the integral (8) and existence of the inverse matrix B−1. The determinant of B,

D2det(B)2=1+t2μxy2μxμy2+t2μx+μy2,
(14)

does not vanish for any dispersion relation Ω(k) and, thus, the problem of wave caustics never appears within this approach as it generally occurs in the lowest orders of the widely used stationary phase and steepest descent methods.

Finally, the reference solution can be written as follows:

η(x,t)=D(k0,t)1/2expΛ(k0,ξ,t)D2/2×cos𝜗0(k0,x,t)+𝜗1(k0,t)+𝜗2(k0,ξ,t).
(15)

We used “mixed” notations introducing “fast” x and “slow” ξ(11) variables in (15) in order to discriminate dependence on the wavetrain carrier parameters (“fast” x) and transformation of the wave packet shape associated with small deviations from the harmonic k0 (“slow” ξ).

In addition to the carrier wave phase ϑ0(9), two phase shifts ϑ1 and ϑ2 appear in (15). The first one ϑ1 shifts phases of all the wavetrain components by the same value,

𝜗1(k0,t)=12arctant(μx+μy)1+t2(μxy2μxμy).
(16)

Generally, ϑ1(k0, t) is decaying as t−1 is dependent on signatures of the second differentials μx, μy, and μxy. For omnidirectional propagation and in the special case of quasi-dispersion (μxy2μxμy=0, see Ref. 12), ϑ1 tends toward the well-known result of the stationary phase method ϑ1 = ±π/4 (e.g., Ref. 9), thus generalizing this result.

The phase shift ϑ2(ξ, k0, t) depends on both the wavetrain dispersion parameters and “slow” coordinate ξ,

𝜗2(ξ,k0,t)=12D2ξ12t3μyμxy2μxμytμx+ξ22t3μxμxy2μxμytμy2ξ1ξ2tμxy1+t2(μxy2μxμy).
(17)

It vanishes for the wavetrain carrier harmonic (ξ = 0) and changes peripheral harmonics when ξ ≠ 0. For deep water waves at t → +, the lines of constant phases of (17) are hyperbolas, and the phase ϑ2(ξ, k, t) is decaying as t−1.

The argument of the exponential function in (15) preserves the Gaussian shape of the wavetrain modulation. The isolines of Λ(k0, ξ, t),

Λ(k0,ξ,t)=ξ121+t2μxy2+μy2+ξ221+t2μxy2+μx22ξ1ξ2t2μxy(μx+μy),
(18)

are ellipses whose parameters (orientation and eccentricity) vary with time. Equation (18) can be written in the following form:

Λ(k0,ξ,t)=t2ξ1μxyξ2μx2+ξ2μxyξ1μy2ξ12ξ22,
(19)

that clearly illustrates effects of the wave packet rotation and stretching from the initial state ξ12+ξ22=const1 to the asymptotic one (t → +) when ξ1μxyξ2μx2+ξ2μxyξ1μy2=const2. It should be stressed that the isotropy of initial distribution and/or of the dispersion relation does not cancel these effects.12 

The factor D affects the wavetrain modulation in two ways in (15): the magnitude of perturbation decays, but its width in a coordinate space grows to preserve the total energy. In the general case, when the dispersion factor does not vanish, i.e., Ωkl2ΩkkΩll0 (μxy2μxμy0), the wavetrain in two-dimensional media decays as t−1. In the special case of quasi-dispersion, Ωkl2ΩkkΩll=0, the decay becomes slower and follows the one-dimensional asymptotics t−1/2.

Thus, the method of this paper reproduces the well-known results of the stationary phase approach. Moreover, it allows for the realization of Lord Kelvin’s idea “to smooth, to smear”13 formally infinite growth at wave caustics in a consistent way. At the same time, it is not a general method of solution for arbitrary initial (boundary) conditions. We leave the issues of completeness and uniqueness of the solutions beyond our analysis. An important point to emphasize is that the method is multi-dimensional. Thus, it essentially extends similar results for one-dimensional waves14 and captures new physical effects like the wave packet rotation.

The general approach presented above is successfully applied here to the problem of ship waves. Formally, it is valid for narrow-banded wavetrains (6), while the ship waves are not. In this regard, additional consideration of wave kinematics is required to substantiate the method for the specific problem. We start this section with an analysis of ship wave kinematics first developed by Lord Kelvin.1 We reproduce and precise some well-known results on phase distributions treated in previous studies.15 Our analysis within the approach of the reference functions provides a simple and effective tool for quantifying the amplitude distributions in the ship wake, the problem that has many important implications including the case of shallow and intermediate water depth.16 In particular, new analytic results of the section explain recent findings for ship wakes at high speeds.17,18

In the reference frame of the ship moving along the x-axis in a negative direction with a constant speed U, the dispersion relation for deep water waves takes the form of the Doppler-shifted frequency,

Ω(k)=g|k|1/2+kU.
(20)

Assuming ship waves to be stationary with reference to the ship (Ω = 0), one immediately has the conditions of synchronism that do not depend explicitly on the dispersion relation (type of waves),

c=Ucosθ,
(21)

where c is the phase speed of deep water waves, c = g/ω with intrinsic wave frequency obeying the dispersion relation ω=g|k|1/2, and

cosθ=k(k2+l2)1/2,
(22)

where θ is the angle between the opposite direction of axis x and the wavevector direction. Formulas (20)–(22) assume real-valued solutions for the stationarity condition Ω = 0, thus implying the far field approximation of the problem under study.

Our approach defines the carrier wave trajectory as ξ = 0 [(6) and (13)], i.e.,

x=Ωkt;y=Ωlt.
(23)

It is equivalent to the geometric optics approach or the stationary phase method which treats the trajectory as one of the quasi-particles traveling at group speed. For the deep water dispersion relation (20) with conditions of stationarity [(21) and (22)], it gives a parametric form of the trajectory,

x=Ut112cos2θ;y=Ut2cosθsinθ.
(24)

Within the problem of stationary ship wake, the time t can be treated as a formal parameter of the stationary phase approach. Indeed, looking for wave phase extrema in (1) with the constraint of wave system stationarity Ω(k) = 0 (20), we obtain the time t as the Lagrangian multiplier and Eq. (21) for the conditional phase extrema.15 Alternatively, within our approach, the carrier harmonic with k and x is fixed by the reference solution shape and does not imply solution of the extremum problem.

The trajectories of (23) are straight lines, and the phase of the carrier harmonic ϑ0 of (9) along these trajectories can be written in a trivial way,

𝜗0=gt2c.
(25)

It leads to the well-known parametric solutions of Kelvin–Lighthill1,5 in terms of phase ϑ0 and angle θ,

x=U2𝜗0gcosθ(2cos2θ);y=U2𝜗0gsinθcos2θ.
(26)

Accordingly, lines of constant phases ϑ1 and ϑ2 [(16) and (17)] can be parameterized by excluding dependence on parameter t.

Turning points of the curves (26) are determined by stationary points, where dx/ = 0, dy/ = 0,

cos2θ1=2/3,
(27)

i.e., for angle θ1 ≈ 35.26°. These turning points form a cone with a spread angle

θK=arctany1/x1=arctansinθ1cosθ12cos2θ1=arcsin(1/3)19.47,
(28)

which is the well-known angle of Kelvin.

One can easily show that the stationary point in coordinates (x, y) at the Kelvin angle θK(28) corresponds to the inflection point of isofrequency (21) at angle θ = θ1(27) in a wavenumber space (k, l), i.e., 2k/l2|θ=θ1=0. The latter condition means merging of two wave systems: diverging and transverse at their demarcation line. The inflection point 2k/l2|θ=θ1=0 corresponds to maximal group velocity when stationary waves can exist and to the maximal opening angle, the Kelvin angle θK, in the coordinate plane (Fig. 2).

Ship wave kinematics is illustrated by two figures. Figure 1 shows wave parameters in a wavevector space, while Fig. 2 describes the phenomenon in the coordinate plane. In Fig. 1(a), two isofrequency curves Ω(k) = 0 [see Eqs. (20) and (21)] are shown for different dimensionless ship speeds and Froude numbers Fr = 0.25; 1.0. The latter can be naturally defined in terms of the bandwidth Δk or perturbation scale Δx,

Fr2=U2gLx=U2gΔx=U2Δkg;  Lx=Δx.
(29)

It should be emphasized that wave kinematics of the carrier wave determined by trajectories (23), i.e., by formally point-wise initial conditions, does not depend on the Froude number (29). The corresponding dependencies appear in dispersive corrections of the wave phases ϑ1 andϑ2 [(16) and (17)] and, what is more important, in expressions for wave amplitudes (15) that essentially account for the finite size of an initial perturbation.

The inflection point of the dispersion relation (20) [black dotted lines in Fig. 1(a)] demarcates two systems of ship waves: transverse (shown in blue) and diverging (red) ones. At any coordinate point inside the Kelvin cone [arctany/x<θK], two solutions of the condition of synchronism (21) exist that correspond to one transverse (θ < θ1, blue) and one diverging wave (θ > θ1, red) with corresponding wavevectors (shown in the same color) and different phase shifts ϑ0, ϑ1, andϑ2 [see general expressions (16) and (17)]. Evidently, wavevectors of each family are the same at any point of a wave trajectory (24), shown by a green line in Fig. 2. These wavevectors (shown by colored lines in Fig. 2) are perpendicular to isophase curves ϑ0 = const of transverse and diverging waves. The corresponding identical direction of group velocities of both waves (αT = αD) for the particular trajectory (green line) is shown by a black arrow in Fig. 2. Note that absolute values of group velocities are different and, hence, phases ϑ0, ϑ1, and ϑ2 of diverging and transverse waves in a point are, generally, different.

Two points of the analysis of Sec. III A should be emphasized. First, all the wave trajectories x(t), y(t) as functions of t are straight lines, and there are no intersections between them (caustics). Second, at every point of these trajectories, two harmonics (transverse and diverging) obeying the synchronism condition (21) coexist. Thus, wave interference at any point of the stationary ship wake is represented by two harmonics only. The stationarity condition (21) wipes out all the satellites of the carrier wave making the wave packet narrow. In this way, the inherently wide-banded problem of ship waves is reduced to a problem of interference of just two narrow-banded wavetrains. Thus, the formulas for the evolution of narrow-banded perturbation can be modified by formal passage to δ-pulse reference wavetrains (3).

In fact, key expressions for wavetrain evolution (15)–(19) contain the bandwidth parameters Δk, Δl explicitly [(6) and (8)],

μx=(Δk)2U34gcos3θ(3sin2θ1),μy=(Δl)2U34gcos3θ(3cos2θ1),
(30)
μxy=3(ΔkΔl)U34gcos4θsinθ,μxy2μxμy=(ΔkΔl)2U6cos6θ8g2.
(31)

These small but still finite quantities (Δk and Δl) should be taken into account as they essentially determine wave amplitudes via the determinant D(14).

The bandwidth parameter Δkx) can be associated with the dimension of the traveling ship and should satisfy, formally, the condition of the narrow-banded perturbation that matches the stationarity condition (21). The latter can be re-written for the wavevector modulo æ =(k2+l2)1/2 and its maximal perturbation Δæ,

Δææ=UcosθgLx=Frcosθ1.
(32)

Simple algebra for determinant D in (14) gives

D2=1+Fr82β2𝜗02cos8θ2+Fr84𝜗02cos8θ(β21)(3cos2θ1)+12,
(33)

where the parameter

β=Δx/Δy
(34)

is a length–beam ratio of the ship.

Bandwidth narrowness is provided by Gaussian distribution (4) that can be re-written taking into account the stationarity condition [(21) and (22)] as follows:

F(k,l)=expFr42(1cosθ)2cos2θ+tan2θβ2cos2θ.
(35)

The representation (35) illustrates an effect of amplitude (energy) partitioning between different wave systems that depends on the angle parameter θ and the wavevector direction (22). Figure 1(b) shows the correspondence between the ranges of wavenumbers (frequencies) and the wave amplitude content of different wave systems. Waves with relatively low transversal wavenumbers (transversal waves) gain energy from the core of distribution F(k, l) (4), while diverging waves with high transversal components of wavevectors are controlled by low energy tails of F(k, l). Thus, at low Froude numbers Fr, transverse waves (relatively small θ) predominate, while diverging waves with large angles θ (1/cos θ ≫ 1) are suppressed.

At the Kelvin cone cos2θ = 2/3 that demarcates diverging and transverse waves, the “weight” of the corresponding harmonics is written as follows (35):

FKelvin(Fr,β)exp0.05+0.75β22Fr4,
(36)

which means an effective suppression of diverging waves at small Froude numbers. For Fr = 0.3 and β = 5 (length–beam ratio), wave amplitude at the edge of the Kelvin cone appears to be less than 1% of the maximal amplitude of the transverse waves at the wake axis. For Fr = 0.5, this fraction becomes about 53% and decays slower inside the Kelvin cone. This can explain the observed growth of the diverging fraction of the ship wake at high Froude numbers and the visible tendency of the ship wake to shrink at high Froude numbers Fr ≫ 1. The latter effect is interpreted as shrinking of the Kelvin cone,19 while in terms of our approach is naturally treated as energy redistribution between transverse and diverging wave fractions inside the permanent Kelvin cone.

The limit of small Froude numbers (32) does not simplify analysis. Moreover, the Froude number Fr enters D in (33) in combination with parameters β and ϑ0 that could be and are indeed large values. Thus, asymptotics

Fr8𝜗02β2+;Fr8𝜗02+,
(37)

make sense and can be taken into consideration as well.

An additional phase shift ϑ1 associated with the wave dispersion (16) depends on the bandwidth Δkx or Fr),

𝜗1(k0,t)=12arctanFr4𝜗0cos4θ(β21)(3cos2θ1)+12+Fr8β2𝜗02cos8θ.
(38)

At small Froude numbers (32) and finite β, ϑ0, it gives the result by Lighthill,5 

lim(Frcosθ)0𝜗1=0.
(39)

For large values of the length–beam ratio β and ϑ0, i.e., at the limit (37), a new characteristic value of angle cos2θ2 = 1/3 appears in addition to cos2θ1 = 2/3. This new angle is implied by symmetry of elliptical initial perturbation in our problem [cf. (28)],

θN=arctany2/x2=arctansinθ2cosθ22cos2θ2=arcsin(233)15.79.
(40)

For the length–beam ratio β2 ≫ 1, the second term in (33) has a minimum near θN ≈ 15.8°, i.e., in the Kelvin cone interior. It can be observed in the ship wave patterns irrespective of the value of Froude number.17,18

At large phases ϑ0, the ship wake amplitude in (15) tends to simple power-like dependence,

AD1/2Fr4cos4θ𝜗01,
(41)

which is also consistent with results of Rabaud and Moisy.17,18

With explicit formulas at hand, one can easily present the role of different mechanisms in stationary ship wake formation. Formally, the approach implies an essential restriction: the bandwidth of an initial perturbation (boundary conditions in the case of stationary ship waves) should be narrow which leads to a formal requirement Fr cos θ ≪ 1 (32). In fact, this limitation is not essential as long as the stationarity condition (21) guarantees the bandwidth narrowness by wiping out satellite harmonics, thus making the wave packet spectrum narrow. In contrast to the general reference solution (15), which depends on both time t and coordinate x, the stationarity condition (21) reduces the number of arguments of the particular solution (33) and, additionally, prunes wave interference in a wavevector space to a superposition of just two harmonics of diverging and transverse waves.

The resulting modulation of the ship wake spectrum by the Fourier image F(k, l) (4) can be formally used in formulas of Sec. II. Figure 3 presents patterns of ship wakes at different Froude numbers Fr = 0.5, 1.0, and 2.0 (from the top to the bottom). The length–beam ratio β = 10 is fixed. The scaling of Rabaud and Moisy18 is used to relate the results of different approaches. The amplitudes are normalized by maximal values for x/Lx > 2. The “near-field” range x/Lx < 2 is left blank.

The left column of Fig. 3 shows diverging waves. The central column shows the transverse wave system, and the right one shows the resulting superposition of transverse and diverging waves. The figure illustrates the mentioned amplitude effects. The diverging waves are quite weak at Fr = 0.5, as expected from formulas (35) and (36). At high Froude numbers (middle and bottom rows in Fig. 3), they become dominating as compared to the transverse wave system. A visible collapse of the ship wake cone at high Froude numbers is associated with the mechanism of energy partitioning between different wave systems, as described and illustrated in Fig. 1(b). This transparent physical treatment of the effect of transformation of the ship wake at high Froude numbers does not require an introduction of long scale cutoff of initial perturbation spectra as proposed by Rabaud and Moisy18 and Darmon et al.19 

In this paper, the reference solution approach for linear water waves is presented. It allows for transparent physical analysis of the classic problem of stationary ship waves. Explicit formulas are obtained for wave variables (free surface elevation in our case) without explicit reference to dynamical equations. All the required information is contained in the fact of wave linearity, i.e., the possibility of superposition of partial solutions and in the specific dispersion relation Ω(k).

The approach allows for separate analysis of amplitude and dispersion effects for independent systems of diverging and transverse waves. A transparent physical treatment of the effect of amplitude partitioning (Fig. 1) is presented and illustrated by wave patterns at different Froude numbers. Snapshots of ship wakes in Fig. 3 obtained with explicit analytical formulas show consistency with the results of previous studies based on conventional asymptotic approaches.17,18 Within the novel analytical approach, the important effect of visible ship wake shrinking at high Froude numbers found its natural physical explanation as a phenomenon of energy partitioning between different families of stationary ship waves (diverging and transverse).

In contrast to the conventional stationary phase and the fastest descent approaches, the two-dimensionality of the problem is not an essential difficulty. Moreover, the resulting formulas allow for easy discrimination between different physical effects associated with wave dispersion. The effect of phase ϑ1(38), which is a straightforward counterpart of the phase shift ±π/4 of the method of the stationary phase and the dependencies of wave parameters on the Froude number and length–beam ratio β have not been detailed here.

The novel approach of reference solutions has good perspectives. In particular, the problem of a ship wake in the presence of shear current20 can be analyzed easily within this approach.

The presented method of reference solutions by Gnevyshev and Badulin12 is associated with a specific class of initial (boundary) conditions and cannot be used in general cases. At the same time, the linearity of the problem essentially extends the approach applicability to qualitative and quantitative analysis as a testbed for numerical simulations.

This research was supported by the Russian Science Foundation (Grant No. 19-72-30028) with the contribution of MIGO GROUP (http://migogroup.ru).

1.
Lord
Kelvin
,
Proc. R. Soc. Edinburgh
25
,
1060
(
1906
).
2.
W.
Thomson
,
Proc. Inst. Mech. Eng.
38
,
409
(
1887
).
3.
T. H.
Havelock
,
Proc. R. Soc. A
81
,
398
(
1908
).
4.
H.
Lamb
,
Hydrodynamics
, 6th ed. (
Dover
,
New York
,
1932
),
738
p.
5.
J.
Lighthill
,
Waves in Fluids
(
Cambridge University Press
,
Cambridge, United Kingdom
,
1978
),
504
p.
6.
W. R.
Hamilton
, “
Researches respecting vibration, connected with the theory of light
,” in
The Mathematical Papers by Sir William Rowan Hamilton
, Dynamics Vol. II (
Cambridge University Press
,
1940
), Chap. XIX, Note Book 52, pp.
451
526
.
7.
W.
Thomson
,
London, Edinburgh, Philos. Mag. J. Sci., Ser. 5
22
(
137
),
353
(
1886
).
8.
M. V.
Fedoryuk
,
Asymptotic: Integrals and Series
(
Nauka
,
Moscow
,
1987
).
9.
M. V.
Fedoryuk
, “
Saddle point method
,” in
Encyclopedia of Mathematics
(
Springer Science+Business Media B.V.; Kluwer Academic Publishers
,
1994
), ISBN: 978-1-55608-010-4.
10.
G.
Stokes
, “
Problem 11 of the Smith’s Prize examination papers (Feb. 2, 1876)
,” in
Mathematical and Physical Papers
(
Cambridge University Press; reprinted by Johnson Reprint Co.
,
New York
,
1966; 1905
), Vol. 5, p.
362
.
11.
V. P.
Maslov
,
Asymptotic Methods and Perturbation Theory
(
Nauka
,
Moscow
,
1988
).
12.
V. G.
Gnevyshev
and
S. I.
Badulin
,
Moscow Univ. Phys. Bull.
72
,
415
(
2017
).
13.
W.
Thomson
, “
On ship waves
,” in
Popular Lectures and Addresses: In Three Volumes
, Nature Series Vol. III (
MacMillan and Co.
,
London; New York
,
1891
), pp.
450
506
.
14.
G. F.
Clauss
and
J.
Bergmann
,
Appl. Ocean Res.
8
,
190
(
1986
).
15.
P. N.
Svirkunov
and
M. V.
Kalashnik
,
Phys.-Usp.
57
,
80
(
2014
).
16.
T.
Soomere
, “
Long ship waves in shallow water bodies
,” in
Applied Wave Mathematics
(
Springer-Verlag, Berlin-Heidelberg
,
2009
), pp.
193
228
.
17.
M.
Rabaud
and
F.
Moisy
,
Phys. Rev. Lett.
110
,
214503
(
2013
).
18.
M.
Rabaud
and
F.
Moisy
,
Ocean Eng.
90
,
34
38
(
2014
).
19.
A.
Darmon
,
M.
Benzaquen
, and
E.
Raphaël
,
J. Fluid Mech.
738
,
R3
(
2014
).
20.
S. Å.
Ellingsen
,
J. Fluid Mech.
742
,
R2
(
2014
).