Laboratory investigations on the resonant feature of `dead water' phenomenon

Interfacial internal wave excitation in the wake of towed ships is studied experimentally in a quasi-two layer fluid. At a critical `resonant' towing velocity, whose value depends on the structure of the vertical density profile, the amplitude of the internal wave train following the ship reaches a maximum, in unison with the development of a drag force acting on the vessel, known in the maritime literature as `dead water'. The amplitudes and wavelengths of the emerging internal waves are evaluated for various ship speeds, ship lengths and stratification profiles. The results are compared to linear two- and three-layer theories of freely propagating waves and lee waves. We find that despite the fact that the observed internal waves can have considerable amplitudes, linear theories can still provide a surprisingly adequate description of subcritical-to-supercritical transition and the associated amplification of internal waves.


Introduction
When a ship is travelling through a strongly stratified water body a certain amount of its kinetic energy is being used up for the excitation of internal waves in its wake, hardly noticeable from the surface, yet perceived as a drag force acting on the vessel. For centuries, the phenomenon has been known to Norwegian seamen as "dödvand" or dead water. In the fjords of the Scandinavian coastline fresh water from slow glacier runoff gently sets on top of the saline seawater without substantial mixing and hence nearly jump-wise vertical density profiles can develop [1]. These circumstances facilitate particularly strong dead water effect associated with large amplitude wave activity along the internal density interface (pycnocline).
In the logbook of the 1893-96 Norwegian Polar Expedition Arctic explorer Fridtjöf Nansen reported experiencing marked dead water drag on board research vessel Fram that reduced the ship's speed to a fifth part. Nansen's original observations were further analyzed by (the then-PhD student) V. W. Ekman, who, in order to understand the phenomenon, has conducted laboratory experiments in a quasi-two-dimensional wavetank, filled up with a two-layer working fluid consisting of saline-and freshwater. In these measurements a scale model of Fram was towed along the surface subjected to either constant or gradually changing force as control parameter, against which the model's velocity U was measured [2]. Ekman found that for smaller towing forces (and smaller U ) dead water drag F dw follows a quadratic scaling F dw ≈ ζU 2 up to a certain threshold, where coefficient ζ is a function of the density profile ρ(z) and the wetted area of the ship.
However, when U reaches a critical value of U ≈ 0.8·c (2) 0 , where c (2) 0 denotes the long-wave velocity of in-arXiv:1901.05931v1 [physics.flu-dyn] 17 Jan 2019 Fig. 1 Hysteresis of ship speed -as observed and discussed by Ekman in his original work -illustrated graphically as a function of changing towing force (decreasing branch: blue, increasing branch: red). The unstable branch, where the dead water phenomenon occurs, is marked with a dashed line. The green curve represents the (inverse of) subcritical relationship F dw ≈ ζU 2 .
terfacial waves on the pycnocline (to be discussed later) the corresponding F dw starts to decrease significantly. Then, in Ekman's "changing force"-type experiments a hysteresis was encountered, as sketched in Fig. 1. When the gradually increasing towing force reached a certain tipping point (the right end of the blue branch in Fig.  1) the speed of the ship suddenly jumped to a higher value. An analogous drop of U could be observed when decreasing the towing force in time (red branch of Fig.  1). The unstable branch (dashed line) is inaccessible to the system for any prescribed towing force. When initiating experiments with a constant force within this hysteretic regime, Ekman found that the ship's velocity U exhibited large fluctuations, comparable to the mean value [2].
If one intends to explore the dynamics in this unstable branch it is therefore beneficial to conduct experiments in which, instead of the applied force, the towing velocity U is prescribed. Such settings are common in the laboratory modeling of lee-wave dynamics (see, e.g. [3,4,5]), where obstacles of various shapes are being towed at the surface (or at the bottom) in a tank of stratified working fluid and the properties of the generated internal waves, hydraulic jumps and 'wave rotors' are evaluated [6,7,8].
In the present experimental work we follow a similar approach to address the dead water phenomenon, by pulling a ship model at a constant speed and -instead of the experienced drag -analyzing the properties of the excited interfacial waves on the pycnocline. Here the 'critical regime' is characterized by wave patterns whose vertical extent is comparable to the height of the upper layer. The fact that wave amplitudes show strong dependence on the wave velocity (that is set by ship speed U ) implies that the waves are of nonlinear nature [9]. Even for velocities where no wave trains are observed, the localized bump following the ship at the pycnocline resembles the solitary wave solutions of the nonlinear Korteweg-de Vries equation [10,11].
A large pool of theoretical, numerical and experimental studies exists discussing the properties of internal waves emerging in the dead water problem addressing nonlinear wave excitation both in the subcritical [12] and supercritical [13] regimes, as well as the applicability of the theoretical findings to observational data, including the historical logs of Fram [14].
Despite of the obvious nonlinearity of the problem, as far as the interfacial wave forms are concerned, the effect of nonlinear corrections is found to be surprisingly minor (see, e.g. Fig. 5 of [13]). In a series of earlier laboratory experiments on topography-induced large amplitude interfacial waves [15] we have also demonstrated that linear three-layer theories, e.g. the one of Fructus & Grue [16] may yield remarkably good fits to experimental data. Therefore, in the present work we focus on the applicability of linear theories to the dead water phenomenon in similarly stratified settings. Our aim here is to explore the U -dependence of the wavelength and amplitude of internal waves in the wake at fixed towing velocities and contrast the results with predictions of the linear theory for lee waves and for freely propagating three-layer interfacial waves. These measurements supplement the 'constant force' experiments of Mercier et al. [17] who have conducted state-of-theart measurements with a presribed towing force (utilizing a falling weight) to propel a ship model in a rather similar set-up.
The paper is organized as follows. Section 2 describes the experimental set-up and the applied data acquisition methods. Our results are presented in Section 3. The paper is then concluded with a brief discussion of the findings in Section 4.

Experimental set-up and measurement methods
The experiments reported here have been carried out in a rectangular laboratory tank made of transparent plexiglass. Its length and width are L = 239 cm w = 8.8 cm, respectively (see Fig. 2). The tank was filled up to level H = 12 cm with density stratified water: the bottom domain contained saline water solution colored by red or blue food dye for the sake of visualization. Following the preparation of this layer, freshwater was poured through a sponge slowly onto the water surface in order to minimize mixing effects and to yield quasi-two layer density profiles, characterized by an approximately 2 cm thick region of steep density increase (see Fig. 3) or 'gradient layer'. The temperature differences within the water body were negligible.
The properties of the 11 different stratification profiles are summarized in Table 1. h (2) 1 and h (2) 2 mark the effective thicknesses of the upper and bottom layers, respectively, using the two-layer approximation (hence the upper index '(2)'), i.e. assigning a jump-wise effective density change to the mid-height of the gradient zone. As a characteristic vertical scale in the system a 'reduced thickness' H r can be introduced as the harmonic mean of h 2 ). Table 1 also lists the values of the average density ρ 2 of the bottom layer and the corresponding twolayer linear interfacial wave velocity c (2) 0 in the longwave limit that reads as where g denotes the gravitational acceleration and ρ 1 ≈ 1 kg/l is the average density in the top (freshwater) layer [18].
For a more precise treatment of the density profiles a three-layer approximation can also be applied [16], in which the top, gradient and bottom layers (indexed with j = 1, 2, 3, respectively) are characterized by their approximate thicknesses h (3) j and density gradients, or equivalently, their buoyancy (or, Brunt-Väisälä) frequencies N j that take the form The values of h (3) j and N j are obtained via piecewise linear regression to a given profile by determining the Fig. 2 The schematics of the set-up. The geometrical parameters of the tank are L = 239 cm (total length, not fully shown), w = 8.8 cm, and H = 12 cm.
intersection points and slopes of the fitted lines. In the three-layer theory there is no such an explicit formula for the long wave velocities c (3) 0 as in the two-layer approximation (1), as will be addressed later. Hence, the c Table 1 are numerical results.
To capture dead water phenomenon we investigated the interfacial internal wave excitation behind towed LEGO TM 'tug boat' models (series 4005 and 4025) [19]. These toy ships have a modular design of a bow and a stern building block and four removable identical intermediate segments. Thus five configurations (S 1 to S 5 ) with different lengths d could be investigated, as listed in Table 2 and shown in Fig. 4. The width of the models was 5.9 cm, comparable to tank width w. The ship model was towed by a (nylon) fishing line spanned horizontally above the water surface by ca. 10 cm, as sketched in Fig. 2, and driven by a DC motor whose voltage (and hence, the towing speed) could be adjusted between the experiment runs. The ship models' draught (i.e. the vertical distance between the waterline and the bottom of the hull) was found to be approximately 1 cm in all cases.
Each experiment was recorded with a HD video camera (at frame rate 50 fps and frame size 720px × 1280px) pointing perpendicularly to the sidewall close to the middle of the tank. To acquire the precise value of ship velocity U and to obtain time series of vertical motion of the interface the video recordings were evaluated by Tracker, an open source correlation-based feature tracking software [20].

Qualitative description of the flow
As a ship model moves along the tank in the studied velocity range U ∈ (1.3; 12.2) cm/s it generates pronounced waves on the internal interface, while the displacement of the free water surface remains negligible, as visible in Fig. 4. An important property of the observed dynamics is that the internal waves are following the ship and propagate at the same velocity as the ship itself. This is visualized in the space-time plots of Fig. 5 for three different constant towing speeds U (see caption). In these diagrams the shading of a point at horizontal position x and time t is given by the sum darkness of the pixel column at x as calculated from the video frame at time t (e.g. the ones shown in in Fig. 4). Since the background of the tank is stationary throughout the videos, the spatial and temporal changes in darkness are attributed to internal waves. The trajectory of the bow of the ship is highlighted with solid black line in each panel.   Table 2 The lengths of the used ship configurations.
The typical wavelength λ and the characteristic amplitude A of the internal waves are set by the towing speed U , the density profile ρ(z), and the ship's length d. The A(U ) dependence is far from monotonous: for each stratification there exists an intermediate velocity U * at which internal waves of the largest amplitude develop. This 'resonant' amplification is demonstrated with the video frames in Fig. 6 for three values of towing velocity U , listed on the panels (stratification #11, ship configuration S2). The snapshots are aligned such that the largest displacements of the density interface (marked with black vertical lines) are beneath each other for better comparability. The ship model's direction of motion is leftward in all images.
At the smallest U (uppermost panel) a wave crest appears at the interface beneath the ship model. For larger U 's a train of waves form (cf. Fig. 5) and the first crest shifts towards the stern of the ship, while the characteristic horizontal size (peak-to-peak wavelength λ) of the interfacial disturbance increases. Coincidentally, its vertical size (or, amplitude A) also increases and reaches a maximum at U * , as captured in the second panel of Fig. 6. In the U > U * regime wavelength λ(U ) continues to increase, but amplitude A(U ) starts to decrease, as seen in the two bottom panels in Fig.  6. In the following subsection we explore the parameter dependence of U * and the associated internal wave dynamics.

Parameter dependence of the resonance
The values of U * were found to range between 2.1 cm/s and 9.1 cm/s in the different experiments. Based on ear- Fig. 4 The five different ship configurations of increasing length d, marked S1 to S5 (downward. cf. Table 2). The bottom layer is visualized using red food dye. The snapshots are organized such that the first wave trough locations behind the ship models are underneath each other.
lier results [21,22] the relevant nondimensional velocity scale of the dead water problem is the internal Froude number i.e. the ratio of the ship speed U and the twolayer long wave velocity c (2) 0 . Plotting the largest observed vertical interface displacement (peak amplitude) A against U/c (2) 0 for each experiment indeed yields a fairly good data collapse, as shown in Fig. 7a with a peak at approximately 0.8. Each data point represents a single experiment run (towing) and the different symbols (and colors) indicate the various ship configurations used, as shown in the legend.
Empirical 'resonance curves' have been fitted to the above data for all five ship configurations separately, yielding the solid curves of Fig. 7a (their coloring matches that of the corresponding data points). The function takes the 'usual' form of Parameter a measures the maximum of the curve, b represents the 'width' of the resonance and C accounts for the position of the peak. The fitted values of b and C exhibit clear monotonous dependence on ship length d as shown in Figs. 7b and c. The b(d) dependence even appears to be consistent with a linear scaling. The error bars represent the regression errors of the given parameters in both panels. It is to be noted that for all ship configurations the data points in Fig. 7a manifestly scatter in a wider range around the fit in the U > U * (supercritical) range than for lower ship speeds. This implies that in this regime some physical parameter of the stratification profiles other than the Froude number U/c (2) 0 becomes relevant, as will be addressed later.

Comparison with linear two-and three-layer theories
In what follows we compare the observed wave (and ship) speeds U and wavenumbers k to the available theoretical predictions of linear two and three layer theories.
Assuming two homogeneous water layers of different densities separated by a sharp interface, the phase velocity reads as where the notations are as introduced in Section 2 [23].
In the long-wave (k → 0) limit the relationship takes the form of equation (1), therefore c (2) (0) ≡ c 0 . Expressing the velocities and wavenumbers in the problem's 'natural' nondimensional units, i.e. U/c (2) 0 and kH r (hereafter referred to as k ), respectively, maps equation (4) to the same invariant curve for all twolayer density profiles. This graph is shown with blue solid line in Fig. 8, alongside the measured data points. The symbol shapes mark different ship configurations (see legend) and the coloring represents nondimensional wave amplitude A , i.e. the maximum vertical displacement A of the interface divided by the parameter a of the fitted resonance curve of the given ship configuration (cf. Fig. 7). k was calculated via measuring peakto-peak wavelengths λ = 2π/k between the second and third wave troughs. Note, that this method yields a considerable smaller value of λ than the typical length of the first wave trough (cf. Fig. 5). Runs where no wave train developed were omitted from this analysis.
As expected, the linear two-layer theory systematically overestimates the wave speeds for larger values of k . Interestingly, however, the best agreement between the data and the theory (in the k ≈ 0.5 domain) is observed at near-resonant wave speeds, i.e. where wave amplitudes A are large (see Fig. 8).
For larger towing speeds (U/c (2) 0 0.8) the observed wavenumbers remain larger than the predictions of the two-layer theory and exhibit a roughly inversely proportional scaling U/c (2) 0 ∝ k −1 (see the gray hyperbolic guide curves in Fig. 8) the implications of which will be discussed in subsection 3.4.
The linear three-layer theory described in [16] has provided a quite accurate description of the observed interfacial velocities and wavenumbers in flow-topography  interaction experiments [15]. This approximation assumes rigid top surface, small-amplitude waves, and three layers with depths h (3) j (j = 1, 2, 3) and piece-wise linear stratification characterized by buoyancy frequencies N j (see eq. 2 and Table 1). The dispersion relation c (3) (k) can then be derived numerically from the implicit equation Table 1 corresponding to the long-wave (k → 0) limit are also derived numerically from the above formulae.
To demonstrate the difference between the predictions of the two-and three-layer theories, the rescaled c (3) (k) curve calculated with the parameters of stratification profile #5 (see Table 1) is added to Fig. 8 in the form of a red curve. Note, that unlike the aforementioned two-layer curve (blue) the three-layer one is not invariant at all in the units used here. For instance, in this exemplary case c holds, but for many other profiles the sign would be reversed (cf. Table 1).
In the correlation diagrams of Fig. 9a and b we plot the theoretical phase velocities of the two-and threelayer models, respectively, against the measured wave speeds U for each observed wavenumber k. In both panels the theoretical long wave velocity (c was used as the unit for nondimensionalization. Symbol shapes denote different ship configurations and the coloring represents rescaled amplitude A as in Fig. 8. As noted before, the two-layer theory systematically overestimates the speeds in the U/c (2) 0 0.8 range (i.e. when U U * ), thus the vast majority of the data points scatter above the y = x line (black) in panel a. The three-layer theory, however, yields a fairly good match with the observations in the same subcritical regime. In the supercritical range, however, both two-and threelayer approximations break down entirely, as indicated by the deviation of the data points from the y = x line, implying that here another physical mechanism, besides  rescaled with the 2-layer long wave velocity, given by eq. (1) -as a function of nondimensional wavenumber k = kH r . The symbol shapes mark different ship configurations (see legend), and the color scale marks non-dimensional amplitude A , rescaled by the fitted maximum values of the corresponding resonance curves (cf. Fig. 3). The blue curve represents the invariant (with respect to the nondimensional units used here) two dimensional velocity-wavenumber relation (4), the red curve denotes an exemplary three-layer relation, obtained for the parameters of experiment #5, based on the implicit formula 5. The grey curves represent two hyperbolae,i.e. 'isofrequency curves' in the corresponding nondimensional time units, y = 0.5/x and y = 1/x, respectively. the ship speed U , becomes relevant in the wavenumber selection.

Lee-wave dynamics
A hyperbola in a velocity-wavenumber dispersion plot marks a constant frequency. For the nondimensional parameters of Fig. 8 the hyperbolic guide curves represent identical frequencies with respect to the stratificationdependent time unit H r ρ 1 /(g(ρ 2 − ρ 1 )). Thus, the fact that the data points appear to follow hyperbolic scaling when U > U * implies that a characteristic frequency associated with the given density stratification ρ(z) determines the observed wavenumbers.
Such physically meaningful 'eigenfrequencies' are the buoyancy frequencies N j of the layers (see Table 1), among which the mid-layer value N 2 is the largest in all cases. Fitting the function α/k to the dimensional U (k) data in the U > U * (supercritical) range yields the empirical frequency parameter α that is plotted against N 2 in Fig. 10. The error bars represent the regression errors and, as before, the different symbols mark various ship configurations. The scattering of the data points indicate a linear relationship α = 0.48(±0.01)N 2 . The result of the fit is shown with a black solid line.
The simplest (linear) theory that describes fixed frequency wave propagation behind an obstacle in a stratified fluid is that of lee waves emerging at the down-  stream sides of mountain (or seamount) ridges [7,8].
If the system is characterized by a single buoyancy frequency N (linear stratification) and the obstacle is moving at velocity U with respect to the medium then U = N cos(φ)/k lee holds [23], where φ is the tilt of the wavenumber vector k lee from the horizontal (|k lee | = k lee ). Thus the maximum achievable k lee becomes k lee = N/U , corresponding to fully horizontal wave propagation.
As mentioned above, buoyancy frequencies N 1 and N 3 of the top and bottom layers are small in the studied density profiles, thus any lee wave activity is expected to be confined to the roughly 2 cm-thick intermediate layer which would then act as a 'waveguide' [15] for the observed oscillation frequencies (i.e. above N 1 and N 3 ) at the interface. Since the layer thickness is an order of magnitude smaller than the typical wavelengths in the U > U * regime, k lee is nearly horizontal here (φ ≈ 0). An important difference from the 'traditional' lee wave setting is that in our configuration the ship model does not reach down to the intermediate layer, thus it is the depression of the interface underneath the ship and not the model itself that acts as a 'moving obstacle' for the flow in the intermediate layer. The linear scaling presented in Fig. 10 appears to be consistent with the proposition that it is indeed lee wave-like dynamics that is observed in the supercritical regime, but -due to the aforementioned peculiarities of the stratification profile -not in its 'classical' form. As an empirical correction we may therefore introduce an effective buoyancy frequency N eff ≈ 0.48N 2 to characterize the stratification, whose physical interpretation will be addressed in the Discussion.
Finally, we investigated the critical wavelength λ * (≡ 2π/k * ) corresponding to the maximum amplitude as a function of ship length d. The results are shown in Fig.  11 for seven different stratification profiles, marked by different colors and symbols. The solid black line marks y = x and the gray guides are linear functions with a slope of 0.5 and a vertical offsets 10 cm and 27 cm. Apparently, the critical wavelength tends to increase with ship length, that roughly follows the empirical formula λ * ≈ 0.5d + f (ρ(z)), implying that the ship configuration also plays a role in the wavenumber selection.

Discussion and conclusions
Inspired by the historical work of Ekman analyzing the 'dead water' phenomenon [2] we have conducted laboratory experiments on wave excitation by a ship model towed at a fixed speed over a salinity stratified water body. We have analyzed the dependence of the interfacial wavenumber k and peak-to-peak amplitude A on the towing speed U , the stratification profile ρ(z) and ship length d.
Due to the fact that in this setting the observed amplitudes are comparable to the characteristic vertical length scale H r of the problem, the excited waveforms can only be explained, if at all, by means of nonlinear wave theories (see, e.g. [12,13]). Yet, we deliberately focused our analysis on linear approximations in order to explore to what extent can these account for the basics of the observed dynamics, most notably the resonancelike amplitude amplification around U/c (2) 0 = 0.8 and the associated transition in terms of the U (k) dependence.
From the findings reported in the previous section it appears that the observable characteristic wavenumber k at a given U is set by the larger one of the corresponding 3-layer wavenumber k (3) (U ) predicted by the linear approximation of Fructus and Grue [16] and the lee wavenumber k lee (U ) derived using 'effective buoyancy frequency' N eff ≈ 0.48N 2 . In other words, among the two competing mechanisms the one yielding shorter waves generates the first trough behind the ship model and, hence, sets the characteristic wavelength in the system. The crossing point of the two dispersion relations where k lee = k (3) holds is encountered around the critical towing speed U * , as sketched in Fig. 12. At this resonant wavenumber the two wave types would be superimposed onto each other resulting in an amplified interfacial wave excitation, that is confirmed by the observations. As a secondary effect, the selection of the resonant wavelength λ * was also found to be influenced by the length d of the ship, as shown in Fig. 11 in unison with the increasing trend in the U * (d) dependence (Fig. 3c).
The reason for the occurrence of the aforementioned factor of 0.48 in the effective buoyancy frequency N eff is unclear. Comparing this value to the internal wave dispersion relation ω = N 2 cos(φ) within the intermediate layer we find that it would imply a wave propagation whose lines of constant phase lie at an angle φ ≈ 61 • to the vertical. However, we were not able to identify any geometrical constraint (e.g. one associated with the depression at the interface behind the ship) that would necessitate the presence of such a limiting angle.
On the other hand, taking the total bottom-to-surface density difference ∆ρ and assuming linear stratification along the full depth H gives a certain 'mean buoyancy frequency' N m = g/ρ 1 (∆ρ/H) which is found to follow the average scaling N m ≈ 0.46(±0.04)N 2 for the profiles listed in Table 1, that is fairly close to N eff . Thus, the explanation may be that for long waves whose wavelength λ exceeds H, the otherwise complex threelayer profile can be simply 'averaged over' in terms of density gradients.
Our results clearly demonstrate the somewhat surprising and unexpected message that linear theories can occasionally be applied for the description of interfacial waves in such velocity and amplitude ranges that other- Fig. 12 The schematic wavenumber velocity diagrams of the lee waves (red curve), the 3-layer model of Fructus and Grue [16] (black), and the actually observable U (k) domain (green) set by the maximum of the two competing U (k) functions at each k. The intersection point U * (k * ) (black filled circle) fairly corresponds to maximum amplitude internal wave excitation.
wise belong to the domain of nonlinear wave dynamics. It is clear, however, that linear models are insufficient to describe more complex features of the observed phenomena, e.g. waveforms, vorticity, etc. The authors also hope to increase awareness in the community about the dead water phenomenon, which -despite being discovered over a century ago -is still an interesting ground for theoretical, numerical, and experimental research.