The dead water problem, in which under certain conditions a vessel advancing in a stratified fluid experiences a considerable increase in resistance respect to the equivalent case without stratification, was studied using computational fluid dynamics (CFD). The advance of a vessel in presence of a density interface (pycnocline) results in the generation of an internal wave that in the most adverse conditions can increase the total resistance coefficient by almost an order of magnitude. This paper analyses the effects of stratification on total and friction resistance, the near field wake, internal and free surface waves, and resistance dynamics. Some of these effects are reported for the first time, as limitations of previous efforts using potential flow are overcome by the use of a viscous, free surface CFD solver. A range of densimetric Froude numbers from subcritical to supercritical are evaluated changing both the ship speed and pycnocline depth, using as platform the research vessel athena. It was found that the presence of the internal wave causes a favorable pressure gradient, acceleration of the flow in the downstream of the hull, resulting in thinning of the boundary layer and increases of the friction resistance coefficient of up to 30%. The total resistance presents an unstable region that results in a hysteretic behavior, though the characteristic time to establish the speed–resistance curve, dominated by the formation of the internal waves, is very long and unlikely to cause problems in modern ship speed controllers.
A surface vessel sailing in stratified seas with a speed close to the celerity of the internal waves results in a phenomenon known as “dead-water,” which causes the boat to experience increased resistance and loss of steering power. The internal wave-making resistance is analogous to the free surface wave-making resistance, and it can result in forces several times higher than those for a vessel sailing in nonstratified waters.
Internal waves are gravity waves affecting a stratified fluid, that it is a fluid with vertically varying density. Internal waves are volumetric in nature, but can be represented in terms of deformation of the interface between two fluids, or if the density varies continuously with depth of the location where vertical density gradient is maximum, called the pycnocline, can generally be used to define the position of the interface. In the upper part of the ocean layers can be caused by large differences in either salt concentration (for instance, due to fresh water from glaciers or littoral streams above the denser ocean salt water) or temperature. A strong and stable density stratification can be obtained after long periods with very calm seas, a rationale for the term dead-water. A ship sailing through this environment produces internal waves with a Kelvin pattern similar to the ones observed in surface waves. The energy needed to form these internal waves results in a large change in resistance, and consequently propeller thrust. If the required thrust to overcome such resistance exceeds the vessel engine capacity, the maximum speed of the vessel can be severely limited.
The dead-water phenomenon was first reported by Nansen , whose ship experienced increased resistance in the Arctic Ocean during his journey to the North Pole. Ekman  conducted an early experimental study by towing a ship model in a tank with stratified conditions, concluding that the model experienced an additional drag due to the internal waves, and that the drag force peaked as the ship speed reached the celerity of the interfacial waves.
Two-layer models of the stratification are frequently used to study the phenomenon, since most physical features are present with this simpler representation of the stratification. Hudimac  used a moving point source to develop a formula to calculate the wave-resistance of a thin ship for both internal and surface waves. He observed that divergent and transverse waves are present for speeds lower than a critical value, but only divergent waves were present for larger speeds. Miloh et al.  used linear theory to compute the motion of a prolate spheroid in a two-layer fluid with finite depth, focusing on wave resistance and internal and free-surface disturbances. They found a sudden rise of internal wave resistance at a densimetric Froude number . Other studies [5,6] reported similar findings using experimental or theoretical techniques. Of particular interest is the work of Mercier et al. , which studied three-layered and linear stratification, showing that the dead-water phenomenon is still present in these more realistic situations, and of Grue , who showed that under subcritical speeds the dead-water resistance coefficient only depended on the densimetric Froude number and that the maximum value depended on the ship draft.
This paper analyzes the dead-water phenomenon using computational fluid dynamics (CFD) for the Research Vessel Athena operating in a two-layered fluid. The study includes for the first time results of effects on friction resistance, not available from analyses using potential flow. The dynamics of the formation of the internal waves are analyzed in the context of the time evolution of the resistance and the dynamics of the vessel, discussing the conditions at which the vessel can overcome the increased resistance and the dead-water speed limit, as well as the presence of a hysteretic behavior in the speed–resistance curve.
Mathematical and Numerical Modeling
Variable density effects have been fully incorporated to an existing naval hydrodynamics solver, maintaining the incompressibility assumptions (i.e., pressure changes do not affect density distribution). This approximation is appropriate for all situations in which density changes are mostly caused by concentration of external agents like temperature and salinity. Density effects have been incorporated to all conservation equations, including the turbulence model. A full description of the mathematical model can be found in Ref.  and references therein, and only a brief summary is presented here.
The calculated level set is used to impose appropriate free surface boundary conditions for the other fields in the single-phase simulations, where only the water field is solved. Details on the implementation and numerical requirements can be found in Ref. .
The mathematical model described has been implemented in the computational fluid dynamics suite REX. REX is a high-order, structured, multiblock URANS solver with delayed detached eddy simulation capabilities, developed primarily for naval hydrodynamics, but that has also been applied successfully to other fields. Pressure and velocity fields are coupled using a modified projection method , designed to improve continuity conservation. REX uses overset technology  to provide connectivity between moving blocks, which allows arbitrary motions of the grids defining the geometry at different hierarchical levels, from complex maneuvers of the vessel, to motions of controlled appendages and the operation of fully discretized propellers. Proportional–integral–derivative controllers are implemented for heading, speed, and standard maneuvers [15,16]. An engine model has been implemented for more realistic simulation of conditions for which thrust output is a limiting constraint for the vessel.
For all the results presented herein, second-order Euler backward difference was used for the temporal terms in all equations. Convection terms used hybrid high-order upwind schemes (second to fourth). Variable number of inner iterations was used per time step, depending on convergence of the solution, to a maximum of four. Standard boundary conditions are imposed for the transport of the density field, employing user-provided profiles for inlets and also as initial condition.
Figure 1 shows a schematic of the simulation conditions. The upper layer has density and depth (the location of the pycnocline), with the lower layer of infinite depth and density . The far-field upstream density profile is imposed as:
which represents an interface that transitions from to in with smooth derivatives and a maximum dimensionless density gradient .
The U.S. Navy research vessel Athena in barehull form and in even keel condition is chosen as the ship geometry. R/V Athena has been subject to extensive research work including simulations under several conditions in model and full scale  and references therein. Main particulars are shown in Table 1.
The overset grid system is shown in Fig. 1. Only the starboard side of the ship was simulated, taking advantage of the centerplane symmetry of the problem. The grid consists of two blocks, the hull and background grids. Three systematically refined grids with grid refinement ratio were generated to evaluate grid convergence. The resulting coarse, medium, and fine grids have 7.7 M, 21.9 M, and 61.5 M grid points, respectively. The simulations required for the grid study were run with different time steps so that the Courant number is the same, resulting in , and for the coarse, medium, and fine grids, respectively. Only two cases were run using all three grids, due to computational cost, while all other cases were run with the medium grid. The total number of points for all runs was between 17 and 22 million points for the medium grid, depending on the position of the density interface and expected internal wave amplitude. The system is distributed for parallel computation to between 69 and 89 processors, with all computations conducted at the Helium and Neon clusters at the University of Iowa. Nonslip boundary conditions are imposed at the hull, while inlet and exit are imposed at the upstream and downstream faces of the background grid, respectively, centerplane symmetry at , and farfield conditions on the other faces. Initial and inlet conditions impose also the density stratification profile in Eq. (7).
Simulations are performed in model scale with the hull fixed in even keel condition, though results are reported in dimensionless form or, when relevant, in full scale using only Froude scaling with no Reynolds number corrections. The Froude, Reynolds, and bulk Richardson numbers are , and , respectively, computed with reference length and speed of (the ship length) and (47 m and 6.25 knots in full scale). Nonstratified simulations set and no stratification in the initial and boundary conditions, but are otherwise identical to simulations with stratification.
Results and Discussion
Ship Moving in a Two-Layer Stratified Fluid.
The resistance induced by the internal wave peaks at , with . Thus, the same Richardson number (or densimetric Froude number) can be achieved varying the depth of the pycnocline ( in dimensionless form, with the draft of the ship) or the ship speed. The effects of these three parameters on resistance are analyzed next.
Simulations were conducted in three modes: (1) varying densimetric Froude number in the range by changing pycnocline depth () and maintaining speed at ( knots), (2) varying densimetric Froude number in the range by changing speed ( to , 1.14–4.61 knots) at a constant pycnocline depth , and (3) maintaining constant densimetric Froude number at by changing speed and pycnocline depth simultaneously in the range . Resistance coefficients for all cases are shown in Fig. 2, which also includes the resistance coefficient for nonstratified flow. Electronic Annex I, which is available under the “Supplemental Materials” tab for this paper on the ASME Digital Collection shows a video illustrating the slow evolution of the system, due to the typical times for develop the internal wave system and reach steady-state.
Two speeds at constant pycnocline depth are used to evaluate grid convergence, resulting in and . Figure 3 shows the time histories of resistance for these conditions. In the case of lower speed, the densimetric Froude number is close to critical where the resistance changes rapidly with (see Fig. 2), and there is a significant difference between coarse and medium grids, but the difference between medium and fine grids is very small. The coarse grid is unable to properly resolve the internal wave, particularly at the stern, but the medium and fine grids predictions are very close, as shown in Fig. 4. At higher speed, the resistance due to the internal wave is less sensitive to , as shown in Fig. 2, and the three grids predict similar results.
The resistance coefficient shown in Fig. 2 offers rich information. The resistance peaks between and , depending on the location of the pycnocline, and decreases rapidly for lower and higher . These results are in agreement with the study of Miloh et al.  for a prolate spheroid, which peaked at . Notice also that the resistance coefficient is about six times higher than for the case without stratification.
For computations maintaining constant pycnocline depth or ship speed, the resistance coefficient is approximately the same at subcritical regardless of whether the pycnocline depth or the speed is varied. This is expected as a subcritical corresponding to strong stratification conditions that dominate over other effects. However, at larger , the resistance coefficient is higher when the same densimetric Froude number is achieved decreasing the pycnocline depth instead of increasing the ship speed.
Figure 5 shows the variation of the resistance coefficient with pycnocline depth for different densimetric Froude numbers. The resistance coefficient peaks for pycnocline depths in the range , conditions where the density interface is touching the bottom of the hull. This indicates strong near-field effects of the interactions between the hull and the internal wave, which will be discussed in Sec. 4.3. Using a strongly nonlinear potential flow method, Grue  reported similar trends for the effects the draft of a hullform acting on an interface at fixed depth.
Of particular interest is the change in friction resistance in stratified flow, which cannot be predicted with potential flow methods. Figure 6 shows that increases in friction resistance are considerable, and as it is the case with the total resistance coefficient, peak at Frh = 0.8–0.9. The friction resistance peaks for pycnocline depths in the range for subcritical , with much smaller influence of interface depth at supercritical . This increase in frictional resistance is a result of thinning of the boundary layer caused by a favorable pressure gradient and local flow acceleration due to the constraint that the density interface presents to the flow. This phenomenon will be further discussed in Sec. 4.4. Note that the viscous component contributes with 16% of the total resistance at the peak for stratified flow, while for nonstratified flow the friction accounts for approximately 80% of the total resistance.
As has been well established by other researchers , the internal waves are characterized mainly by the densimetric Froude number and dimensionless obstacle height , as shown in Figs. 7 and 8, and display some similar characteristics to free surface waves. The highest amplitude in Fig. 7, which shows the internal wave for , occurs at , coincident with the peak in resistance coefficient in Fig. 2. The occurrence of the peak at subcritical has been previously reported [4,19], and corresponds to dispersion effects due to the finite depth of the top layer. Notice also that a distinct wake in the internal wave is observed, another characteristic that cannot be captured with potential flow methods. At subcritical densimetric Froude numbers ( and ), the presence of transverse waves is clear. behaves as a near critical densimetric Froude number, slightly lower than the expected . This behavior has been reported also in the experiments of Lacaze et al. , who attribute the difference to the lack of excitation of faster, long waves than the ship length that would be able to travel upstream.
Figure 7 also shows the predicted wavelength and wake half angle for Kelvin waves at the lowest subcritical condition , and the angle of the divergent waves at the supercritical condition , as discussed in Tulin et al. . While in ships the wave pattern typically exhibits bow and stern Kelvin waves at low Froude numbers, the internal wave interacts with the hull in a different manner, creating a more complex wave system, but still maintaining the 19.47 deg half angle in the line joining the first line of wave crests, and a wavelength close to the predicted for Kelvin waves.
At the supercritical condition, , the transversal waves disappear and a Mach cone formed by divergent waves is observed , with an angle close to the 52 deg predicted in Ref.  for these conditions. The effect of depth on the internal waves is shown in Fig. 8 for . The wave field far from the ship is fairly independent of the depth of the pycnocline, but effects are strong near the body. Most noticeably, immediately downstream of the stern the internal wave forms a pattern similar to a “rooster tail” for surface waves, shown in the inserts and most noticeably for . This pattern is typical of dry conditions in transom stern ships like R/V Athena operating at high Froude numbers . At higher , a rooster tail is still visible but separates farther downstream of the stern and produces two peaks and a trough.
At higher pycnocline depths, the recovery point of the rooster tail occurs further downstream, causing the crest locus of the stern divergent wave to move downstream as well. In the wake of the ship, a stationary pattern is observed, similar to the one reported by Lacaze et al. . Higher densimetric Froude numbers result in narrower internal wave wakes, see Fig. 8, as is also the case in surface waves, with higher pycnocline depths having the same effect.
The free surface waves, created by the interaction of the air/water interface with the advancing ship, are of very small wavelength and height due to the low Froude numbers ( at which the ship operates when presenting dead water effects (the grid density used in this paper is too coarse to resolve these very small free surface waves), and in nonstratified flow, the free surface is essentially flat. In stratified flow, the free surface mirrors the internal wave with opposite sign, and peaks in the free surface appear in troughs in the internal wave and vice versa. This occurs because the flow is accelerated in restricted areas (peaks in the internal wave), reducing the pressure and sucking in the free surface. The opposite occurs in troughs in the internal wave. Figure 9 shows that the amplitude of the free surface elevation is approximately 2.5% that of the internal wave for and .
The presence of the pycnocline imposes interfacial forces that strongly affect the flow. Not considering the small mixing caused by turbulence, the sharp interface effectively acts as a slip boundary condition for the flows in both sides of it. The large-amplitude internal wave occurring for cases close to peak resistance influence the flow field significantly. Figure 10 shows the dimensionless cross () and axial () velocities on the interface for the case with and . The presence of the internal wave causes considerable inflow in the stern of the ship, where the crest of the internal wave is, and restricts the passage area, resulting in an acceleration of the flow and higher axial velocities than for the nonstratified case. This can be seen in Fig. 11 as a low pressure, high-velocity region below the ship in the downstream half of the hull. This higher velocity and resulting favorable pressure gradient result in a much thinner boundary layer than for the nonstratified flow case (for which a scalar with the same initial distribution as density is passively transported by the flow to produce an equivalent iso-surface), with the consequent increase in friction resistance shown in Fig. 6. Similar effects were also reported for a submarine near a density interface by Martin et al. . It is interesting to note that the presence of stratification flattens the interface in the wake, both due to the effect of restitutive buoyancy forces and inhibition of turbulent mixing and stirring, resulting in more consistent vortex detachment from the transom stern than for nonstratified flow. The higher axial velocity in the stratified case also results in higher vortex detachment frequency. The vortex detachment dynamics plays an important role in the behavior of the wake of the ship.
Wave and Resistance Dynamics.
The resistance reported in Sec. 4.2 was obtained after running long enough to achieve either a steady-state or a statistically stationary value. In reality, a simulation where the ship speed is imposed suddenly to a system in rest involves an initially very large force as the geometry overcomes the virtual mass acceleration. There is a transient as the boundary layer develops and the internal wave system evolves. Figure 12 shows the evolution of the dimensionless resistance at different for , highlighting the transient nature of the process. At low densimetric Froude number, the resistance coefficient slowly decreases as approaches steady-state. The opposite occurs at , where resistance grows slowly over time. Higher or intermediate densimetric Froude numbers attain final resistance coefficients faster, but still take 10–20 ship lengths to converge. After a relatively fast transient resulting from the formation of the boundary layer (occurring in 2–3 ship lengths), the transient behavior is dominated by the formation of the internal waves, which is a slow process.
To evaluate the impact of the internal wave dynamic response to the operation of a surface craft, consider a constant thrust propulsion system, with dimensionless thrust . Figure 13 shows the resistance/velocity evolution for runs at constant dimensionless thrust , starting from rest with a pycnocline depth of . The open symbols in Fig. 13 represent the steady-state resistance for each speed, which exhibits a local maximum of at knots, corresponding to in Fig. 12. If the internal wave formed quickly, the ship should never be able to overcome the peak of resistance, since the imposed thrust is lower. However, the internal wave takes some time to form, and as it does, the resistance increases. Thus, a slow acceleration will result in higher resistance and a consequent added difficulty to overcome the dead-water effect.
In Fig. 13, three accelerations at constant thrust are shown, resulting from three different inertia values: the actual inertia of the ship, 3.3 times the inertia of the ship and 10 times the inertia of the ship. Only the inertia in the direction of motion has been altered by increasing the mass of the ship, causing slower acceleration; the corresponding vertical force (the weight) has not been modified, and if the ship was free in the vertical direction (which is not), it would sink to a new position of equilibrium with a deeper draft. The lowest acceleration case follows closely the steady-state resistance points and is unable to overcome the dead water resistance. In this case, the resistance is lower than for faster accelerations at lower speeds, consistent with Fig. 12 that shows that the resistance decreases as the internal wave develops at low densimetric Froude numbers, and the resistance is higher at higher speeds, where approaches . Faster accelerations result in higher resistance at lower speeds and lower resistance at higher speeds, overcoming the peak in resistance at , transitioning to the unstable region of the resistance curve and accelerating until the resistance again equals the thrust at a higher densimetric Froude number. The internal wave fields, when the ship reaches , are shown in Fig. 13, with the fully developed internal wave field corresponding to an infinitely slow acceleration (equivalent to steady-state at every ship speed), the slowest acceleration considered, for which the internal wave field is similar to the steady-state near the vessel but has not been able to fully radiate to the whole domain, and a high acceleration, where the internal wave field is barely starting to form.
The presence of the unstable negative slope resistance region in Fig. 13 causes hysteresis in the speed/resistance curve when operating at constant thrust, as well as problems maintaining speed between 1.5 and 2 knots. In this region, a higher speed results in lower resistance, which causes acceleration of the ship unless the thrust is reduced. Conversely, lower speeds result in increased resistance, and deceleration of the ship unless the thrust is increased.
Figure 14 shows the resistance/speed curves at resulting from the step increases and decreases in thrust shown in Fig. 15. In Fig. 14, the empty symbols indicate the steady-state resistance at each speed, and the small full symbols indicate the points indicated in Fig. 15. Starting from , the dimensionless thrust is increased to then , and , waiting for steady-state to be reached before the next increase, simulating a slow increase in thrust. For thrust below the peak in steady-state resistance at , the resistance/speed curve follows well the steady-state curve (points a, b, c, and d in Fig. 15). When the thrust is increased beyond the peak in steady-state resistance in the transition from points d to e, the ship experiences initially a higher resistance than the steady-state resistance curve, since the internal wave needs some time to evolve to the new instantaneous speeds and was initially in the configuration offering the highest wave resistance. The time evolution of the resistance and internal wave are presented in the video in Electronic Annex II, which is available under the “Supplemental Materials” tab for this paper on the ASME Digital Collection. The transition from points d to e runs through the unstable region of the resistance curve, resulting in a considerable jump in velocity from 1.5 to 2.5 knots.
A series of step decreases in thrust, staring from steady-state at , and going through then, , and , results in similar behavior compared to increasing thrust. Prior to reaching the local minimum steady-state resistance at approximately 2 knots, the resistance evolves following the steady-state resistance curve. When the thrust is decreased from points d′ to e′, the resistance cannot follow the unstable region of the curve and jumps from a velocity of 2–1.45 knots. The time evolution of the resistance and internal wave are presented in the video in Electronic Annex III, which is available under the “Supplemental Materials” tab for this paper on the ASME Digital Collection. Note that since the transient resistance is higher than the final steady-state resistance, the ship speed undershoots its final value in point e′. Note in Fig. 15 that for each thrust increasing the internal wave development and achieving steady-state takes several minutes in full scale, including the transition from points d to e. This jump through the unstable region of the resistance curve is so slow that a speed controller can easily maintain velocity even though it would require constant action.
The dead water problem was studied using CFD. Simulations were performed for the research vessel athena for a range of subcritical and supercritical densimetric Froude numbers, with results of total resistance coefficient similar to those already reported by other researchers using potential flow techniques, with peak resistance occurring for around 0.9. It was shown that changes in pycnocline depth affect significantly the resistance coefficient, with resistance peaking for . The viscous resistance is also affected by the presence of the internal wave, as the favorable pressure gradient that causes the increase in pressure resistance in the most adverse condition also accelerates the flow, resulting in thinning of the boundary layer and a 30% increase in the frictional resistance coefficient. The localized peak in resistance caused by the internal wave results in an unstable region which causes a hysteretic behavior in the resistance/speed curve when operating at constant thrust. The characteristic times of the resistance/speed are very slow, in the order of minutes, dominated by the time required to form the internal waves, thus making maintaining speed with an active controller possible.
Computations were performed at The University of Iowa HPC clusters Helium and Neon.
Office of Naval Research (Grant No. N000141310351).
- Aw =
ship wetted area,
celerity of the longest internal waves, ;
densimetric Froude number;
undisturbed distance between free surface and pycnocline,
turbulent kinetic energy,
- L0 =
reference length, chosen as the ship length in this work,
dimensionless piezometric pressure;
effective Reynolds number;
bulk Richardson number;
source term in turbulence equation
- U0 =
- Us =
- v =
- vturb =
turbulent eddy viscosity,
- λ =
wave elevation respect to undisturbed pycnocline (for internal waves) or surface (for surface waves)
level set function
specific turbulence dissipation rate,