Hostname: page-component-586b7cd67f-t8hqh Total loading time: 0 Render date: 2024-11-24T01:33:58.321Z Has data issue: false hasContentIssue false

The influence of static versus dynamic pressure distribution strategies for modelling nonlinear waves generated by ships

Published online by Cambridge University Press:  22 October 2024

Jinyu Yao
Affiliation:
State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University, Shanghai 200240, PR China
Harry B. Bingham
Affiliation:
Department of Civil & Mechanical Engineering, Technical University of Denmark, 2800 Lyngby, Denmark
Xinshu Zhang*
Affiliation:
State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University, Shanghai 200240, PR China
*
Email address for correspondence: [email protected]

Abstract

A moving static pressure distribution is commonly used to simulate a travelling ship. However, the ship movement changes the fluid velocity around the hull, inducing pressures on the hull surface that are no longer equal to the static pressure. Therefore, we introduce a dynamic pressure correction strategy, which can accurately simulate the impact of the ship movement on the hull-surface pressure and preserve the desired hull shape under both stationary and transient conditions. The strategy is applied to a high-order spectral model and used to investigate ship-induced waves and wave resistance over a both flat and variable topography. We explore various parameters in our study, including the average water depth to ship draft ratio ($h_0/d$), the channel width to ship width ratio ($W/B$), the Froude number ($Fr_0=U/\sqrt {gh_0}$) and variations in bathymetric slope. Compared with experiments on a flat bottom, the numerical results with dynamic correction show better accuracy in the simulation of ship-induced waves and wave resistance than those obtained using a static pressure distribution. The correlation coefficient for wake waves between the numerical and experimental results is improved by approximately 0.25 with the dynamic correction strategy. The amplitude and wavelength of ship-induced mini-tsunamis over a variable topography are found to be reduced when employing a dynamic correction compared with a static pressure distribution, and this effect becomes more pronounced with higher Froude number. The static pressure approach is shown to allow large deformations of the desired hull shape and changes in ship volume which are responsible for the different wave patterns from the two approaches.

Type
JFM Papers
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
© The Author(s), 2024. Published by Cambridge University Press.

1. Introduction

Ships are pivotal in maritime trade, serving as a fundamental means of transportation. At high speeds, ships produce wake waves and upstream waves characterized by substantial amplitudes. Large waves generated by ships can pose threats to the safety of human and marine structures. Additionally, they can contribute to the erosion of seabeds and coastlines. Accurately simulating moving ships and ship-induced waves is an effective way to evaluate the influence of ship-induced waves on marine structures and environments, and prevent harm induced by these waves.

As a ship travels at a constant speed, the resulting wake waves comprise diverging waves and transverse waves, collectively forming a V-shaped pattern. Based on the linear dispersion relation of gravity waves in deep water, Kelvin (Reference Kelvin1887) found that the wake angle induced by a single point source is around $19.47^\circ$. However, the wake angle induced by a ship travelling in deep water is observed to be smaller than Kelvin's classical angle (Dias Reference Dias2014; Wu et al. Reference Wu, He, Liang and Noblesse2019; Pethiyagoda et al. Reference Pethiyagoda, Moroney, Lustri and McCue2021). This phenomenon was demonstrated by Noblesse et al. (Reference Noblesse, He, Zhu, Liang, Zhang, Zhu and Yang2014) using a wave-interference analysis based on a far-field wave pattern that results from the superposition of two Kelvin wakes originating at the bow and stern of a ship. Additionally, Zhu et al. (Reference Zhu, He, Zhang, Wu, Wan, Zhu and Noblesse2015) studied the longitudinal interference effects on the wake angle in finite water depth and found that these interference effects are especially significant in shallow water. They also identified that the effect of water depth is significant, depending on the length-based Froude number and the ratio of water depth to ship length. It should be noted that wave patterns generated by a ship also depend on the depth-based Froude number $Fr = U/\sqrt {gh}$, where $U$ is the ship speed, $g$ is the gravitational acceleration and $h$ is the water depth. In shallow water, the wake angle increases as $Fr$ approaches unity. Once $Fr$ exceeds unity, the wake angle decreases monotonically with further increases in the depth-based Froude number, due to changing wave speed and independent of interference (Pethiyagoda, McCue & Moroney Reference Pethiyagoda, McCue and Moroney2015).

When a ship travels with a near-critical ($Fr\sim 1$) or supercritical ($Fr>1$) speed, upstream nonlinear solitons are generated. Various nonlinear models, such as the Korteweg-de Vries and Boussinesq equations, have been used to study the generation mechanism and stability of these solitons. Using the forced Korteweg-de Vries equations, Wu (Reference Wu1987) found that a well-balanced interplay between wave dispersion and nonlinearity contributes to the generation of solitons. When sailing in a channel, the ratio of the cross-sectional areas of the ship and the channel is a critical parameter influencing the generation of solitons (Gourlay Reference Gourlay2001). In a channel with restricted width, solitons with a mean free-surface level above the still water surface are periodically generated at the bow, which induces a periodic variation of the wave resistance. Solitons with a straight crestline propagate ahead of the bow at a supercritical speed and begin to break at $Fr=1.2$ (Ertekin, Webster & Wehausen Reference Ertekin, Webster and Wehausen1984Reference Ertekin, Webster and Wehausen1986). To accurately simulate wave breaking at the crests or steep fronts of solitons, the shallow-water equation scheme and the eddy viscosity-based scheme were adopted by Shi et al. (Reference Shi, Malej, Smith and Kirby2018) in a fully nonlinear Boussinesq model. A new modified generalized Boussinesq equation was proposed by Li & Sclavounos (Reference Li and Sclavounos2002) to study solitons in an unbounded shallow-water region. Solitons with a parabolic crestline do not break at $Fr=1.2$, but instead form several parabolic water humps that are blocked ahead of the ship.

The mini-tsunami is a new type of upstream wave observed in a harbour at Flaskebekk, Oslofjorden, Norway (Grue Reference Grue2017). Mini-tsunamis are generated when a ship travels from relatively deep water to a shallow-water region at a subcritical speed. The difference between the water depths in the deep and shallow regions, denoted by $\Delta {h}$, is comparable to the average water depth $h_0$, i.e. $\Delta {h}/h_0\sim {1}$. The Froude number based on the average water depth $Fr_0=U/\sqrt {gh_0}$ ranges from 0.4 to 0.7. The generation mechanism of mini-tsunamis is different from that of nonlinear solitons and can be attributed to the interaction between a travelling ship and a topographic variation along the ship track (Grue Reference Grue2020). Mini-tsunamis observed at Flaskebekk are typically 0.5–1 km long, with a period ranging between 30–60 s. These upstream waves, with a maximum wave height of approximately 1.4 m, extend across a 2–3 km wide fjord and can result in significant coastal erosion. To investigate the wave characteristics of mini-tsunamis under different Froude numbers and channel widths, Grue (Reference Grue2017) proposed a linear model in which a travelling ship is modelled by a moving static pressure distribution. The crestline of mini-tsunamis is a curve in a wide channel but a straight line in a narrow channel. As the Froude number increases, the wave height increases exponentially, while the wavelength decreases. Several short waves, which are induced by the gradients of the ship pressure, propagate behind the mini-tsunamis (Grue, Pedersen & Saetra Reference Grue, Pedersen and Saetra2022). The mean free-surface level of these short waves is above the still water level, which can be attributed to second-order free-surface nonlinearities (Yao, Bingham & Zhang Reference Yao, Bingham and Zhang2023).

Pressure integration along the hull surface is a direct method used to compute wave resistance which plays a crucial role in ship optimization design and reducing energy consumption. The pressure along the surface of a moving ship is different from the static pressure on a stationary ship. This difference results in the sinkage and trim of the ship, which leads to an increase in the wave resistance (Ma et al. Reference Ma, Zhang, Huang, Yang, Gu, Li and Noblesse2017). Pressure integration is very convenient to use in methods where an advancing ship is modelled by a moving pressure distribution (Huang & Wong Reference Huang and Wong1970; Doctors Reference Doctors1972Reference Doctors1993; Yeung et al. Reference Yeung, Wan, Banumurthy, Ham and Lew2008). However, pressure integration alone makes it difficult to distinguish the contributions of upstream waves and wake waves to the total resistance. The wave-cut analysis method is an effective approach for computing the resistance induced by the wake waves, under the assumption that a ship sailing at constant speed generates a stationary wave pattern (Wehausen Reference Wehausen1973). In a study by Janson & Spinney (Reference Janson and Spinney2004), comparisons between the longitudinal wave-cut (LWC) and transverse wave-cut (TWC) techniques revealed that the TWC can be truncated to any length as long as it encloses the entire Kelvin wedge. In contrast, the LWC should be extended as far downstream as possible. As a result, the computational domain needed for a TWC analysis is smaller compared with that required for the LWC. Upstream waves, like periodic solitons, are transient in both the global and ship-fixed reference frames. Hence, far-field methods based on either momentum or energy conservation can be employed to compute the resistance induced by these waves (Amini-Afshar & Bingham Reference Amini-Afshar and Bingham2021; Amini-Afshar, Bingham & Kashiwagi Reference Amini-Afshar, Bingham and Kashiwagi2023). A complete and corrected form of the far-field formulation with body-surface integrals was presented by Kashiwagi (Reference Kashiwagi2022) and shown to be equivalent to the time average of the Lagally theorem.

In numerous past studies concerning ship-induced waves, researchers have represented a moving ship by employing a static pressure distribution that can be readily determined using the formula $p=\rho {g}\eta _s$ where $\eta _s$ is the desired hull shape (Ertekin et al. Reference Ertekin, Webster and Wehausen1986; Wu Reference Wu1987; Tanimoto, Kobayashi & Ca Reference Tanimoto, Kobayashi and Ca2000; Grue Reference Grue2017; Shi et al. Reference Shi, Malej, Smith and Kirby2018). When the ship is stationary, the free-surface shape in the ship region corresponds to the desired hull shape. However, when the ship moves or travels over a depth variation, the pressure along the hull surface is no longer equal to the static pressure, so the ship hull deforms and the volume changes, as can be seen in figures 1, 2, 5 and 6 of Ertekin et al. (Reference Ertekin, Webster and Wehausen1986) for example. This deformation can be expected to influence the accuracy of the simulation of ship-induced waves and may lead to inaccurate assessment of the impact of ship-induced waves on marine structures and the marine environment. This serves as the primary motivation for our current investigation. Here, we propose a dynamic correction strategy for the moving pressure distribution. The aim is to ensure that the free-surface shape in the ship region aligns with the desired hull shape regardless of the ship's motion. To accomplish this, we employ a high-order spectral (HOS) model to simulate the generation and evolution of ship-induced waves. Both pressure integration and a combination of wave-cut analysis and far-field momentum conservation are adopted to compute the wave resistance. These methods agree well and are validated by comparison with experimental measurements. Calculations using the static and dynamic pressure strategies are compared to investigate the influence of a changing hull geometry on the generated waves. This analysis encompasses various factors such as wave amplitude, wavelength and wave resistance. Specifically, the study focuses on assessing the influence of dynamic correction on solitons and wake waves induced by a ship traversing a flat bottom, as well as mini-tsunamis generated by a ship navigating a depth change. These analyses are conducted under varying conditions, including different Froude numbers based on the average depth $Fr_0$, the ratios of channel width to ship width $W/B$, the ratios of average depth to ship draft $h_0/T$ and the slopes of the depth variation $\gamma$. The pressure along the hull surface of an advancing ship is no longer equal to the static pressure due to the ship movement. This dynamic effect on the pressure can be accurately captured by the dynamic correction strategy, which preserves the desired hull shape and volume. Changes in hull geometry and volume with the static pressure strategy are quantified and shown to be responsible for increasingly large differences in wave patterns as the Froude number increases.

The remainder of this paper is organized as follows. The dynamic correction strategy in the HOS model is introduced in § 2. The set-up of the numerical calculation and a detailed description of the test cases are presented in § 3. In § 4, the validation of the numerical method is demonstrated by comparison with the experiments. The effects of dynamic correction on ship-induced waves and resistance over a flat bottom are studied in § 5. The impact of dynamic correction on mini-tsunamis is investigated in § 6. The fluctuations in wave resistance as the ship traverses changes in water depth, as well as the nonlinear impacts on resistance are analysed in § 7. Finally, conclusions are provided in § 8.

2. Dynamic correction strategy in a HOS model

2.1. The HOS model

The HOS model is based on the potential flow formalism, where the fluid is assumed to be homogeneous, inviscid and incompressible, and the flow is assumed to be irrotational. The horizontal plane is defined by axis $\boldsymbol {x}=(x,y)$, and $z$ is the vertical axis pointing upwards. The flow is described by a velocity potential $\phi$, which satisfies the Laplace equation,

(2.1)\begin{equation} \Delta{\phi}=0. \end{equation}

The velocity potential on the free surface is defined by $\phi ^S(x,y,t)=\phi (x,y,\eta (x,y,t),t)$, such that the nonlinear free-surface boundary conditions are

(2.2)$$\begin{gather} {\eta_t}=(1+\boldsymbol{\nabla}\eta\boldsymbol{\cdot}\boldsymbol{\nabla}\eta)\phi_z-\boldsymbol{\nabla}\phi^S\boldsymbol{\cdot}\boldsymbol{\nabla}\eta, \quad {\rm on} \ z=\eta(x,y,t), \end{gather}$$
(2.3)$$\begin{gather}{\phi^S_t}=-g\eta+\frac{1}{2}(1+\boldsymbol{\nabla}\eta\boldsymbol{\cdot}\boldsymbol{\nabla}\eta)\phi_z^2-\frac{1}{2}\boldsymbol{\nabla}\phi^S\boldsymbol{\cdot}\boldsymbol{\nabla}\phi^S-\frac{p}{\rho}, \quad {\rm on} \ z=\eta(x,y,t), \end{gather}$$

where $\eta$ is the free-surface elevation, $g$ and $\rho$ are the gravitational acceleration and water density, respectively. Here $\boldsymbol {\nabla }=({\partial }/{\partial {x}},{\partial }/{\partial {y}})$. The moving pressure distribution $p$ is applied on the free surface, and it is introduced in § 2.2 in detail. Periodic boundary conditions are used in the horizontal plane.

The bottom boundary condition for a variable topography is as follows:

(2.4)\begin{equation} \boldsymbol{\nabla}{\phi}\boldsymbol{\cdot}\boldsymbol{\nabla}{\beta}-\phi_z=0, \quad {\rm on}\ z=-h_0+\beta(x,y), \end{equation}

where the total water depth is written as the sum of the average water depth $h_0$ and a depth variation $\beta$. A flat topography is adopted when $\beta =0$.

The detailed solution process for the velocity potential $\phi$ and the free surface elevation $\eta$ in HOS incorporating variable topography was introduced in Gouin, Ducrozet & Ferrant (Reference Gouin, Ducrozet and Ferrant2017). Here, a brief summary of the solution process is provided.

The velocity potential $\phi$ can be expressed as a power series expansion using components $\phi ^{(m)}$, and is truncated at order $M$, i.e. $\phi =\sum _{m=1}^M\phi ^{(m)}$ (see Dommermuth & Yue Reference Dommermuth and Yue1987; West et al. Reference West, Brueckner, Janda, Milder and Milton1987). Here $M$ is the order of free-surface nonlinearity. To account for the variable water depth, the velocity potentials at each order of $m$ are written as

(2.5)\begin{equation} \phi^{(m)}=\phi_{h_0}^{(m)}+\phi_{\beta}^{(m)}, \end{equation}

where $\phi _{h_0}^{(m)}$ is derived at the average depth $h_0$ and satisfies the boundary condition for a flat bottom; while $\phi _{\beta }^{(m)}$ generated by the variable topography can be further expanded with respect to $\xi =\beta /h_0$ under the assumption that $\xi \ll {1}$, and the expansion is truncated at order $M_b$ ($M_b$ is the order of bottom nonlinearity), i.e. $\phi _{\beta }^{(m)}=\sum _{l=1}^{M_b}\phi _{\beta }^{(m,l)}$. Here $\phi _{\beta }^{(m)}$ satisfies the boundary condition (2.4) for a variable topography and a Dirichlet condition on the still water surface,

(2.6)\begin{equation} \phi_{\beta}^{(m)}=0, \quad {\rm on} \ z=0. \end{equation}

The two potential terms $\phi _{h_0}^{(m)}$ and $\phi _{\beta }^{(m)}$ can be expressed by an eigenfunction expansion of free modes which satisfy the boundary conditions,

(2.7)$$\begin{gather} \phi_{h_0}^{(m)}=\sum_{p}\sum_{q}A_{pq}^{(m)}(t)\frac{\cosh\left[k_{pq}(z+h_0)\right]}{\cosh(k_{pq}h_0)}{\rm e}^{{\rm{i}}({k_{xp}x+k_{yq}y})}, \end{gather}$$
(2.8)$$\begin{gather}\phi_{\beta}^{(m)}=\sum_{p}\sum_{q}B_{pq}^{(m)}(t)\frac{\sinh(k_{pq}z)}{\cosh(k_{pq}h_0)}{\rm e}^{{\rm{i}}({k_{xp}x+k_{yq}y})}, \end{gather}$$

where $k_{xp}=2{\rm \pi} {p}/L_x$, $k_{yq}=2{\rm \pi} {q}/L_y$ and $k_{pq}=\sqrt {k_{xp}^2+k_{yq}^2}$. Here $L_x$ and $L_y$ are the lengths of computational domain in the $x$ and $y$ directions, respectively. The modal amplitude $A_{pq}^{(m)}(t)$ can be solved based on (2.2), (2.3) and (2.4) with $\beta =0$ by a fast-Fourier-transform method. By a Taylor expansion of the bottom boundary condition (2.4) around $z=-h_0$ truncated at order $M_b$, the modal amplitude $B_{pq}^{(m)}(t)$ can be computed as a function of $A_{pq}^{(m)}(t)$ (see Gouin et al. Reference Gouin, Ducrozet and Ferrant2017).

After the vertical velocity on the free surface $\phi _z(x,y,z=\eta )$ is computed based on (2.7) and (2.8), $\phi$ and $\eta$ can be time stepped based on (2.2) and (2.3). A fourth-order Runge–Kutta integrator is adopted for the time stepping. The sawtooth instabilities are mitigated using a Savitzky–Golay (SG) filter. The set-up of the SG filter is outlined in Appendix A.

The in-house HOS code, originally developed for a flat bottom based on the method by Dommermuth & Yue (Reference Dommermuth and Yue1987), has been optimized for computational efficiency by transferring input data to graphics processing units (GPU) and employing GPU-accelerated computing techniques. This code has been successfully applied to investigate phenomena such as freak waves in crossing sea states and four-wave resonant interactions in finite water depth (Liu et al. Reference Liu, Waseda, Yao and Zhang2022a; Liu, Waseda & Zhang Reference Liu, Waseda and Zhang2022b). Additionally, building upon the approach by Gouin et al. (Reference Gouin, Ducrozet and Ferrant2017), the in-house HOS code was adapted for a variable bottom. This updated version has been effectively utilized to study the impacts of free-surface and bottom nonlinearities on mini-tsunamis (Yao et al. Reference Yao, Bingham and Zhang2023).

2.2. Dynamic correction strategy

The travelling ship is modelled by a moving pressure distribution $p(x,y,t)$, which is shown in the last term of (2.3). A static pressure distribution can be computed directly from the desired hull shape $\eta _s(x,y,t)$,

(2.9)\begin{equation} p(x,y,t)=\rho{g}\eta_s(x,y,t)={\rho{g}df(x,t)q(y,t)}, \end{equation}

with

(2.10)\begin{gather} f(x,t)=\left\{ \begin{array}{@{}ll@{}} \cos^2\left[\dfrac{{\rm \pi}\left(\left|x-x^*(t)\right|-\frac{1}{2}\alpha_s{L}\right)}{(1-\alpha_s){L}}\right] & \frac{1}{2}\alpha_s{L}\leq{\left|x-x^*(t)\right|}\leq\frac{1}{2}{L},\\ 1 & \left|x-x^*(t)\right|\leq{\frac{1}{2}\alpha_s{L}}, \end{array} \right. \end{gather}
(2.11)\begin{gather} q(y,t)=\left\{ \begin{array}{@{}ll@{}} \cos^2\left[\dfrac{{\rm \pi}\left(\left|y-y^*(t)\right|-\frac{1}{2}\beta_s{B}\right)}{(1-\beta_s){B}}\right] & {{\frac{1}{2}\beta_s{B}\leq{\left|y-y^*(t)\right|}\leq\frac{1}{2}{B}}},\\ 1 & {{\left|y-y^*(t)\right|\leq{\frac{1}{2}\beta_s{B}}}}, \end{array} \right. \end{gather}

where $L$, $B$ and $d$ are the ship length, width and draft, respectively. The centre of the water plane of the ship is located at $(x^*,y^*)$, which varies as a function of time $t$. The two coefficients $\alpha _s$ and $\beta _s$ determine the shape of the wetted hull surface in the $x$ and $y$ directions, respectively (Ertekin et al. Reference Ertekin, Webster and Wehausen1986; Wu Reference Wu1987; Shi et al. Reference Shi, Malej, Smith and Kirby2018). If the ship length, width, draft and the static displacement $V_0$ are determined, the values of the two coefficients $\alpha _s$ and $\beta _s$ can be calculated based on $V_{0}=\iint {p}/(\rho {g})\,{\rm d}\kern0.06em x\,{\rm d} y$. The ship region is a rectangular region defined by the ship's length and width. The ship pressure is zero outside the ship region.

To include the dynamic effects of a moving ship, a dynamic correction for the pressure distribution is derived as follows. By introducing a local coordinate moving with the ship velocity $U$,

(2.12)\begin{equation} \left. \begin{gathered} x'=x-Ut,\\ y'=y,\\ z'=z, \end{gathered} \right\} \end{equation}

the dynamic boundary condition on the free surface (2.3) is written as

(2.13)\begin{equation} {\phi^S_t}=-g\eta+\frac{1}{2}(1+\boldsymbol{\nabla}\eta\boldsymbol{\cdot}\boldsymbol{\nabla}\eta)\phi_z^2-\frac{1}{2}\boldsymbol{\nabla}\phi^S\boldsymbol{\cdot}\boldsymbol{\nabla}\phi^S+U(t)\phi^S_x-\frac{p}{\rho}, \quad {\rm on} \ z=\eta(x,y,t).\end{equation}

The pressure distribution with dynamic correction, which is used to model a travelling ship, is set to

(2.14)\begin{align} p&=-\rho{g}\eta_s-\frac{\rho}{2}\left[\boldsymbol{\nabla}\phi^S\boldsymbol{\cdot}\boldsymbol{\nabla}\phi^S-\phi_z^2(1+\boldsymbol{\nabla}\eta\boldsymbol{\cdot}\boldsymbol{\nabla}\eta)\right]+{\rho}U(t){\phi^S_x +\rho{g}(\eta-\eta_s)} \end{align}
(2.15)\begin{align} &=-2\rho{g}\eta_s+\rho{g}\eta-\frac{\rho}{2}\left[\boldsymbol{\nabla}\phi^S\boldsymbol{\cdot}\boldsymbol{\nabla}\phi^S-\phi_z^2(1+\boldsymbol{\nabla}\eta\boldsymbol{\cdot}\boldsymbol{\nabla}\eta)\right]+{\rho}U(t){\phi^S_x}. \end{align}

The first part of the pressure in (2.14) simply cancels the dynamic pressure in (2.13), to ensure that $\phi ^S_t=0$ as long as $\eta =\eta _s$. The last term is a penalty term which is added to help enforce the condition that $\eta =\eta _s$, and it tends to improve the stability of the method. Further discussion of this penalty term is provided in Appendix B. We note that this is very similar to the scheme originally proposed by Lindberg et al. (Reference Lindberg, Glimberg, Bingham, Engsig-Karup and Schjeldahl2013).

For clarity, the static pressure distribution in (2.9) is referred to as ‘static pressure’, whereas the ship pressure with a dynamic correction in (2.15) is referred to as ‘dynamic pressure’ in subsequent discussions.

To compute the resistance of the ship, the pressure integration method, wave-cut analysis and far-field methods are adopted. These three methods are described in Appendix C.

3. Set-up for the numerical simulations

3.1. Test cases for solitons above a flat bottom

To investigate the properties of ship-induced solitons over a flat topography, Ertekin et al. (Reference Ertekin, Webster and Wehausen1984) conducted a series of experiments in the Ship Model Towing Tank of the University of California, Berkeley. The dimensions of the tank are $L_x=61$ m by $W=2.44$ m. Some of these experimental cases are used in the present study to investigate the effects of a dynamic correction on ship-induced waves over a flat bottom. A Series 60 ship with a block coefficient $C_B=0.8$ travels along the centreline of the channel. The ship length, width and draft are $L=1.524$ m, $B=0.234$ m and $d=0.075$ m, respectively. The water depth over a flat topography is denoted by $h$. Two different water depths, $h=0.125$ m and 0.15 m are considered. The water depth to ship draft ratios are $h/d=1.67$ and 2, respectively. The Froude number $Fr=U/\sqrt {gh}$ ranges from 0.7 to 1.1. A narrow channel and a wide channel, whose widths satisfy $W/B=5.21$ and 20.85, respectively, are adopted. Because the ratio of the width of the experimental tank to the ship width is $W/B=10.43$, a false wall is placed along the centreline of the experimental tank to ensure the channel width satisfies $W/B=5.21$ in test cases with a narrow channel. In the test cases with a wide channel, a half-ship model is towed along the tank wall so that the channel width is equivalent to $W/B=20.85$. Sketches of the set-up for the experiments for both the narrow and wide channels are shown in figure 1. Wave gauges (WGs) are placed in the channel to record the upstream solitons. The specific locations of these gauges are marked in figure 1.

Figure 1. Sketch of the experiments on solitons over a flat topography in the (a) narrow and (b) wide channels (Ertekin et al. Reference Ertekin, Webster and Wehausen1984). The ship moves along the centreline of the channel with a speed $U{(t)}$. The blue lines indicate the ship. Here $L$ and $B$ denote the length and width of the ship, respectively. Red circles with crosses at the centre represent the WGs.

In the numerical simulation, the ship, modelled by a pressure distribution, starts impulsively and then moves steadily along the $x$-axis with a constant speed. The set-up of the numerical simulation is illustrated in figure 2. In the numerical simulation, the channel width $W$ can be assigned any value, eliminating the need to simulate the false wall and the half-ship used in the experiment. The width of the computational domain is set to match the channel width, specifically $W/B=5.21$ with the number of nodes $N_y=128$ for the narrow channel and $W/B=20.85$ with $N_y=512$ for the wide channel. Because periodic horizontal boundary conditions are used in the HOS model and the ship travels along the centreline of the channel, the two boundaries along the $x$-axis can be regarded as the two channel walls. A detailed illustration about this can be found in Appendix D. The length of the computational domain should be sufficient to prevent waves propagating along the $x$-axis from exiting one boundary and re-entering through the other. The length of the computational domain in a non-dimensional form is $L_x/L=80.63$, with the number of nodes $N_x=8192$. The nonlinear order of the free surface is set to $M=3$, meaning that interactions up to third-order nonlinearities are included. To analyse the soliton characteristics and compare the numerical results with the experimental data, four numerical wave gauges (NWGs) are placed in the computational domain; NWGs 1, 2 and 3 are positioned identically to WGs 1, 2 and 3 in the experiments for both the narrow and wide channels; NWG4 is positioned at the same position as WG4 in the experiments with the narrow channel.

Figure 2. Sketch of numerical simulation for the experiments conducted by Ertekin et al. (Reference Ertekin, Webster and Wehausen1984). The ship moves along the centreline of the channel at a speed $U{(t)}$. The blue lines indicate the ship. Here $L$ and $B$ denote the length and width of the ship, respectively. The four red circles with crosses at the centre indicate the positions of the four NWGs.

When the ship sails over a flat topography with a near-critical or supercritical speed ($Fr\geq 1$), both wake waves and upstream solitons are generated. The wake waves form a stationary wave pattern in the ship-fixed frame, and the resistance induced by the wake waves, denoted as $R_{ww}$, can be computed by the wave-cut analysis. As shown in figure 3, a TWC plane moving with the ship velocity $U$ is positioned aft of the stern. The distance between the wave-cut plane and the stern is set at three times the ship length. However, upstream solitons are periodically generated, and the resistance induced by upstream solitons $R_{uw}$ is computed by the far-field method. A far-field control plane, also moving at the velocity of the ship, is established ahead of the bow. The distance between this plane and the bow is configured to be five times the length of the ship.

Figure 3. The positions of the wave-cut plane and the far-field plane used for calculation of wave resistance. Here $L$ is the ship length, and $U$ is the ship velocity.

3.2. Test cases for mini-tsunamis from a variable bottom

Mini-tsunamis are generated by a ship sailing from deep to shallow water. As shown in figure 4(a), the ship travels along the centreline of the channel which is restricted by the two sidewalls parallel to the sailing track of the ship. The coordinates of the centre of the water plane of the ship are $(x^*,y^*)$, where $x^*$ varies with time, but $y^*$ is constant. The ship initially accelerates from a state of rest and then proceeds to sail at a constant speed. The profile of the ship speed variation is defined by

(3.1)\begin{equation} U(t)= \left\{ \begin{array}{@{}ll} U_0\sin(t/T_a), & {{t/T_a\leq{\dfrac{\rm \pi}{2}}}},\\ U_0, & {{t/T_a>\dfrac{\rm \pi}{2}}}, \end{array} \right. \end{equation}

where $T_a$ is a coefficient related to the acceleration duration, and $U_0$ is the stable travelling velocity after the acceleration phase.

Figure 4. Sketch of the test cases on mini-tsunamis over a variable topography. (a) Plan view: the ship moves along the centreline of the channel. The red square represents the centre of the water plane, whose coordinates are $(x^*,y^*)$. The ship accelerates in the relatively deep region, and moves at a constant velocity through a transition in water depth towards the shallower region. The regions between the two blue dashed lines are the transition zone for the water depth. Here $x_a$ and $x_b$ are the centre positions of the two depth changes. (b) Side view: $d$ is the ship draft. The initial position of the ship is denoted by $x_s$. Here $h_1$ and $h_2$ are the depths in the deep and shallow regions, respectively.

The bathymetry changes solely along the channel direction and remains uniform across the channel width. The bathymetric variation is shown in figure 4(b) and is described by $-h_0+\beta (x)$, where $h_0$ is the average water depth and

(3.2)\begin{equation} \beta(x)=-\tfrac{1}{2}\Delta{h}+\tfrac{1}{2}\Delta{h}\left[\tanh\left(\gamma\left(x-x_a\right)\right)-\tanh\left(\gamma(x-x_b)\right)\right], \end{equation}

where $\gamma$ determines the slope of the depth variation, $x_a$ and $x_b$ denote the centre positions of the two depth changes, respectively. The water depths in the deep and shallow regions are denoted by $h_1$ and $h_2$, respectively, with the depth difference given by $\Delta {h}=h_1-h_2$.

In the numerical simulations involving a variable bottom, the ratio of the ship length to ship width is $L/B=6$, and the ratio of ship width to ship draft is $B/d=6$. The two coefficients governing the shape of the wetted hull surface in (2.10) and (2.11) are $\alpha _s=0.5$ and $\beta _s=0.25$, respectively. In the relatively deep region, the water depth is $h_1/d=11.64$, while in the relatively shallow region, the water depth is $h_2/d=4.36$. Therefore, the average water depth is $h_0/d=8$. The size of the computational domain is $(L_x, L_y)=(300h_0,3h_0)$, and the number of nodes is $(N_x, N_y)=(8192,512)$. Because the width of the computational domain is equal to the channel width, the ratio of the channel width to ship width is $W/B=4$. The centre positions of the two depth changes in (3.2) are $x_a=90h_0$ and $x_b=210h_0$, respectively. The initial position of the ship is $x_s=30h_0$, and the coefficient controlling the acceleration duration is $T_a/\sqrt {h_0/g}=90$. The ship finishes the acceleration before reaching the centre position of the depth change $x_a$, and proceeds through the depth change with a constant speed. Wave characteristics and wave resistance are studied across multiple scenarios, wherein the Froude number based on the average water depth $Fr_0=U/\sqrt {gh_0}$, and the slope of the depth variation $\gamma$ are varied. All of the test cases featuring variable topography are outlined in table 1.

Table 1. Main parameters used in the variable bathymetry test cases. Here $Fr_0$ denotes the Froude number based on the average depth and $\gamma$ controls the slope of the depth variation.

It should be noted that if the order of bottom nonlinearity $M_b=0$, a flat bottom corresponding to the average water depth $h_0$ is used. If $M_b\geq 1$, the bathymetric variation $\beta (x,y)$ is considered, and the corresponding water depth is expressed by $-h_0+\beta (x,y)$. The increase of $M_b$ does not affect the variable bathymetry but decreases the truncation error of the velocity potential, i.e. the numerical results have an exponential convergence rate with respect to $M_b$ (Gouin et al. Reference Gouin, Ducrozet and Ferrant2017; Yao et al. Reference Yao, Bingham and Zhang2023). In all the test cases with a variable bottom in the present study, the amplitude and wavelength of mini-tsunamis increase with $M_b$. As $M_b$ increases from 2 to 3, the increase in the wave amplitude is less than 6.3 %. The increase is less than 2 % as $M_b$ increases from 3 to 5 and is less than 0.2 % as $M_b$ increases from 5 to 6. The increase in the wavelength is less than 0.3 % as $M_b$ increases from 2 to 5. Therefore, within the HOS model, the order of bottom nonlinearity is established as $M_b=5$, which is enough to consider the impact of bottom variation on waves. The numerical results also have an exponential convergence rate with respect to the order of free-surface nonlinearity $M$. In all the test cases, the increase in the wave amplitude can reach 20 % as $M$ increases from 2 to 3, but is less than 3.9 % as $M$ increases from 3 to 4. The increase in the wavelength is less than 5.8 % as $M$ increases from 2 to 3, and is less than 0.75 % as $M$ increases from 3 to 4. Therefore, $M=3$ is chosen within the HOS model.

4. Validation of the numerical model

4.1. Validation of the wave resistance calculation

To estimate the wave resistance of an air-cushion vehicle (ACV), Doctors (Reference Doctors1993) introduced a solution based on potential-flow theory and compared their results with model test experiments. The ACV is represented by a moving static pressure distribution,

(4.1)\begin{align} p_{ACV}&= \frac{p_0}{4}\times\left\{ \tanh\left[\alpha_a\left( \frac{x-x^*}{L}+\frac{1}{2} \right)\right]-\tanh\left[\alpha_a\left( \frac{x-x^*}{L}-\frac{1}{2} \right)\right] \right\}\nonumber\\ &\quad \times\left\{ \tanh\left[\beta_a\left( \frac{y-y^*}{L}+\frac{1}{2} \right)\right]-\tanh\left[\beta_a\left( \frac{y-y^*}{L}-\frac{1}{2} \right)\right] \right\}, \end{align}

where $\alpha _a=24$, $\beta _a=24$ and $p_0=0.0127\rho {g}L$. The ship has $L/B=3/2$ and travels on a flat bottom. The water depth satisfies $h/L=3.58$.

The wave resistance $R_w$ of the ACV is computed at various Froude numbers based on the ship length $Fr_L=U/\sqrt {gL}$ within the range 0.4–1.14. When the ship moves at constant speed, only wake waves are generated, no upstream waves are present, so the wave-cut analysis method is used to compute the wave resistance $R_w$. A TWC plane moving with the ship velocity is positioned after the stern at a distance of $3L$. As shown in figure 5, the wave resistance coefficients $C_w=\rho {g}R_w/p_0^2B$ computed by the pressure integration and wave-cut analysis are compared with the experimental data. The results from the two numerical methods exhibit strong agreement with each other, and they both align with the experimental results.

Figure 5. Wave resistance coefficient $C_w=\rho {g}R_w/p_0^2B$ for an ACV as a function of Froude number $Fr_L=U/\sqrt {gL}$. The ACV has $L/B=3/2$ and travels in open water with $h/L=3.58$. The red circles indicate pressure integration results while the blue crosses show the wave-cut analysis results. The lines are splines through the point values. The black triangles show the experimental data in Doctors (Reference Doctors1993).

4.2. Validation of the dynamic correction strategy

Several experimental cases in Ertekin et al. (Reference Ertekin, Webster and Wehausen1984) are chosen to validate the dynamic correction strategy in the HOS model. The experimental set-up and numerical model configuration have been detailed in § 3.1. Figure 6 shows the data from WG4 in different test cases where the water depth is $h/d=1.67$ and the Froude numbers are 0.8, 1.0 and 1.1, respectively. The outcomes from the NWG4, represented by the blue solid lines, are derived using the dynamic pressure distribution. In each plot, the peak of the first main wave is aligned to ${t'/\sqrt {h/g}}=0$. Behind the main wave, several short waves propagate with a mean free-surface level above zero. The amplitude of the upstream solitons and the mean free-surface level rise with the Froude number. The data from the NWG4 correspond reasonably well to the experimental data, although there are some phase differences observed in the short waves behind the main wave, particularly when the Froude number $Fr$ is large. The phase differences could be attributed to the idealization of the hull shape. The bow and stern shapes described by (2.10) and (2.11) are not identical to the actual bow and stern of the ship model used in the experiments.

Figure 6. Comparison of the wave gauge records computed by the HOS ($M=3$) with the dynamic correction strategy and those from experiments. The water depth satisfies $h/d=1.67$. The Froude numbers $Fr=U/\sqrt {gh}$ are (a) $Fr=0.8$, (b) $Fr=1.0$ and (c) $Fr=1.1$. The blue solid lines are computed results and the black circles are the experimental data from Ertekin et al. (Reference Ertekin, Webster and Wehausen1984).

As shown in figure 7, the amplitude of the first upstream wave generated by the dynamic pressure distribution in the water depth $h/d=2$ is compared with the experimental records. The Froude number varies from 0.7 to 1.1. The numerical results agree well with the experimental records, providing confidence that the dynamic correction strategy in the HOS model accurately captures the underlying physics, despite the geometric simplification of the hull form.

Figure 7. The variation of the non-dimensional amplitude $A/h$ of the main upstream wave with the Froude number. The water depth is $h/d=2$. The black circles represent the experimental records from Ertekin et al. (Reference Ertekin, Webster and Wehausen1984). The blue squares represent the numerical results computed by HOS ($M=3$) with the dynamic correction strategy.

5. Dynamic versus static pressure strategies on a flat bottom

5.1. Wave patterns

To examine the impact of dynamic correction on ship-induced waves over a flat topography, the upstream waves and wake waves generated by the static and dynamic pressure distributions are compared. The ship travels along the centreline of a wide channel ($W/B=20.85$). The water depth is $h/d=1.67$. Two test cases with different Froude numbers, $Fr=1$ and 1.1, are used. The detailed set-up can be found in § 3.1.

Once the ship starts to move with $Fr=1$, the surface elevation along the centreline of the ship region is shown in figure 8. Since the ship is stationary at the initial time $t=0$, the hull shapes modelled by the dynamic and static pressure distributions are identical, and they correspond to the desired hull shape. After the ship moves impulsively, the ship bottom modelled by the static pressure shows notable oscillations, with the mean level of the ship bottom exceeding the initial ship draft. The volume of the ship modelled by the static pressure is increased by approximately 28 % during the initial acceleration of the ship, as shown in figure 9. On the other hand, the ship bottom modelled by the dynamic pressure stays flat and the ship volume is very nearly constant with time. The dynamic correction strategy accurately preserves the desired hull shape both while accelerating and when moving at constant speed.

Figure 8. The surface elevation $\eta$ along the centreline of the ship region at the initial stage of ship movement. The initial position of the ship is at $(x-x_s)/L=0$. Here $d$ denotes the ship draft. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively.

Figure 9. Variation of the ship volume $V$ at the initial stage of ship movement. The initial position of the ship is at $(x-x_s)/L=0$. Here $V_0$ is the original ship volume.

As shown in figure 10, the data of the NWG3 computed by the two types of moving pressure distribution are compared. The WG3 is placed near the channel wall, and its location is indicated in figure 2. To access the accuracy of the waves induced by the two types of moving pressure distribution, the experimental data of WG3 are also shown in figure 10. Significant depressions are observed at $195\leq {t/\sqrt {h/g}}\leq 221$ in figure 10(a) and $186\leq {t/\sqrt {h/g}}\leq 204$ in figure 10(b), respectively. The depression emerges along the ship body, and propagates towards the channel wall. The lengths of the depressions generated by the dynamic pressure align closely with those observed in the experiments, albeit being shorter than those induced by the static pressure. The amplitude of depression in the numerical results is smaller than that observed in the experiments. This discrepancy can be attributed to the idealization of the hull shape. Due the difference in the length of the depression, the phase of the wake waves generated by the static pressure is very different from that in the experiments. However, the length of the depression and the phase of the wake waves induced by the dynamic pressure are much closer to the experimental observations. The amplitudes of the wake waves generated by the static pressure are larger than those induced by the dynamic pressure and those in the experiments, particularly for the waves between $266\leq {t/\sqrt {h/g}}\leq 283$ in figure 10(a). The upstream solitons induced by the two types of pressure distribution and those in the experiments are comparable.

Figure 10. Time series of free surface elevation $\eta$ recorded by WG3. The whole gauge records comprise the upstream waves, the depression and the wake waves. The black dash–dotted line represents the experimental data from Ertekin et al. (Reference Ertekin, Webster and Wehausen1984). The blue and red dashed lines denote the results computed by the dynamic pressure and static pressure, respectively. Here (a) $Fr=1.0$, the depression is between $195\leq {t/\sqrt {h/g}}\leq 221$; (b) $Fr=1.1$, the depression is between $186\leq {t/\sqrt {h/g}}\leq 204$.

To quantify the disparities between the numerical results and the experimental records, a correlation coefficient $\epsilon$ is employed. The correlation coefficient is defined as

(5.1)\begin{equation} \epsilon[\eta_{num}(t),\eta_{exp}(t)]=\frac{\int_t \left[\left(\eta_{num}(t)-\overline{\eta_{num}}\right)\left(\eta_{exp}(t)-\overline{\eta_{exp}}\right)\right] {\rm d}t} {\left\{\left[ \int_t \left(\eta_{num}(t)-\overline{\eta_{num}}\right)^2 \,{\rm d}t\right] \left[\int_t \left(\eta_{exp}(t)-\overline{\eta_{exp}}\right)^2 {\rm d}t \right]\right\}^{1/2}}, \end{equation}

where $\eta _{num}(t)$ is the time series of free surface elevation computed by the numerical methods, $\eta _{exp}(t)$ is that recorded in the experiments. The time-averaged values of $\eta _{num}(t)$ and $\eta _{exp}(t)$ are denoted by $\overline {\eta _{num}}$ and $\overline {\eta _{exp}}$, respectively. The complete gauge records shown in figure 10 encompass the upstream waves, the depression and the wake waves. In the two test cases, the depressions are between $195\leq {t/\sqrt {h/g}}\leq 221$ and $186\leq {t/\sqrt {h/g}}\leq 204$, respectively; the upstream waves are between $169\leq {t/\sqrt {h/g}}<195$ and $160\leq t/ \sqrt {h/g}<186$, respectively; the wake waves are between $221<{t/\sqrt {h/g}}\leq 292$ and $204<{t/\sqrt {h/g}}\leq 292$, respectively. The values of correlation coefficient $\epsilon$ for upstream waves, wake waves and the whole gauge records are shown in table 2. The correlation coefficients of upstream solitary waves exceed 0.94 for both methods, indicating that the upstream solitary waves generated by both methods closely correspond to the experimental observations. However, the wake waves are more sensitive to the dynamic correction. Due to the phase differences and larger wave amplitudes, the correlation coefficient of the wakes waves for the static pressure is lower than 0.4, whereas for the dynamic pressure it exceeds 0.6. The dynamic correction strategy enhances the correlation coefficient for wake waves by approximately 0.25. Moreover, the correlation coefficient of the entire gauge records for the dynamic pressure is also higher than that for the static pressure. Therefore, the dynamic correction strategy improves the accuracy in simulating the ship-induced waves.

Table 2. Correlation coefficients $\epsilon$ between the numerical and experimental results shown in figure 10. The correlation coefficients of upstream waves, wake waves and the whole gauge records for the dynamic and static pressure distributions are calculated, respectively.

5.2. Wave resistance

Periodic upstream solitons and wake waves are generated at near-critical and supercritical conditions ($Fr\geq 1$). The free surface along the centreline when solitons are generated at the bow is shown in figure 11, where the ship travels along the centreline of a narrow channel with a width of $W/B=5.21$, in relative water depth of $h/d=2$, at $Fr=1.0$. The blue and red lines show the dynamic and the static pressure results, respectively. The region between $11.55\leq (x-x_s)/L\leq 12.55$ is the ship region where the ship modelled by the dynamic pressure distribution has a flat bottom, while oscillations are shown on the bottom of the ship modelled by the static pressure distribution. Therefore, with the dynamic correction strategy, the described hull shape is preserved while the ship is moving. At around $(x-x_s)/L=13.56$, the amplitudes of upstream solitons induced by the dynamic and static pressure distributions are comparable. However, differences in both wave amplitude and phase are shown in the wake waves exhibited at $(x-x_s)/L\leq 11.55$. The wake waves induced by the dynamic pressure have a smaller amplitudes than those generated by the static pressure.

Figure 11. The free surface elevation $\eta$ along the centreline of the channel with a width of $W/B=5.21$ at $t=133\sqrt {h/g}$. The water depth satisfies $h/d=2$.

The resistance induced by the wake waves, denoted as $R_{ww}$, and the resistance induced by the periodic upstream solitons, denoted as $R_{uw}$, are computed by the wave-cut analysis and the far-field method, respectively. The positions of the TWC plane and the far-field control plane are given in § 3.1. The total wave resistance can be calculated from $R_w=R_{ww}+R_{uw}$. Correspondingly, the resistance coefficients are expressed by $C_w=\rho {g}R_w/p_0^2B$, $C_{ww}=\rho {g}R_{ww}/p_0^2B$ and $C_{uw}=\rho {g}R_{uw}/p_0^2B$, respectively. The symbol $p_0=\rho {g}d$ where $d$ is the ship draft shown in (2.9).

Because of the periodic generation of upstream solitons, the variation of resistance computed on the far-field control plane over time is also periodic, as shown in figure 12. Here $C'_{uw} (t)$ is the instantaneous resistance coefficient generated by the upstream solitons. The period $T_s$ of the solitons is $44\sqrt {h/g}$, and the discrepancy in resistance over different periods is relatively small. The time-averaged value of $C'_{uw} (t)$ over one period is used as the resistance coefficient induced by the upstream waves, denoted as $C_{uw}$.

Figure 12. The variation of the instantaneous resistance coefficient generated by the upstream solitons $C'_{uw} (t)$ with time. Here $T_s$ is the period of the solitons. The channel width and the water depth satisfy $W/B=5.21$ and $h/d=2$, respectively.

The fluctuation of wave resistances induced by the dynamic and static pressure distributions over one soliton period is compared in figure 13, and the resistance coefficients obtained by the different methods are summarized in table 3. The wake waves form a stable wave pattern in the ship-fixed reference frame, such that the resistance corresponding to wake waves remains relatively constant over time. Due to the larger amplitudes induced by the static pressure distribution, the resistance coefficient $C_{ww}$ is larger than that for the dynamic pressure distribution. Additionally, because the upstream waves generated by the two types of pressure are comparable, the values of $C_{uw}$ are also nearly identical. Therefore, the total wave resistance induced by the dynamic pressure distribution is smaller than that for the static pressure distribution. To validate our results, the total resistance coefficient $C_w$ computed by the pressure integration (PI) method and from experimental (EXP) records are also presented in figure 13 and in table 3. The periodic variation of the instantaneous resistance coefficient $C'_{uw} (t)$ implies that the instantaneous coefficient of the total resistance ${C'_{uw} (t) + C_{ww}}$ also exhibits periodic behaviour. The total resistance coefficient calculated by $C_w=\overline {C'_{uw} (t)+C_{ww}}$ corresponds well with that obtained by the pressure integration method. However, the total wave resistance induced by the dynamic pressure distribution aligns well with the experimental results. This comparison illustrates that a more accurate wave resistance can be obtained through the dynamic correction strategy.

Figure 13. Variation of the wave resistance coefficients over one soliton period $T_s$. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively. The blue solid line represents the resistance coefficient from the wake waves $C_{ww}$. The red solid line represents the coefficient of the instantaneous resistance induced by the upstream waves $C'_{uw} (t)$, and its time-averaged value in one period $C_{uw}=\overline {C'_{uw} (t)}$ is denoted by the red dashed line. The coefficient of the instantaneous total wave resistance computed by $C'_{uw} (t)+C_{ww}$ is denoted by the black solid line, its time-averaged value in one period $C_{w}=\overline {C'_{uw} (t)+C_{ww}}$ is denoted by the black dashed line. The magenta dashed line and the green dash–dotted line represent the total wave resistance coefficients $C_w$ obtained by the pressure integration (PI) method and the experiments (EXP), respectively.

Table 3. Comparison of the resistance coefficients computed by the numerical models and those in the experiments. Here $C_{ww}$ and $C_{uw}$ are the coefficients of resistance induced by the wake waves and upstream waves, respectively. Here $C_w$ is the total wave resistance coefficient. Calculation methods are indicated in parentheses. Here ‘WA’, ‘FM’, ‘PI’ and ‘EXP’ denote ‘wave-cut analysis’, ‘far-field method’, ‘pressure integration’ and ‘experiment’, respectively.

6. Dynamic versus static pressure strategies on a variable topography

6.1. Generation of mini-tsunamis

Figure 14 shows an example of the mini-tsunami generation process using the conditions outlined in Test 3 from table 1. The ship, modelled by a dynamic pressure distribution, travels at $Fr_0=0.5$. Due to the narrowness of the channel, i.e. the ratio of channel width to ship width $W/B=4$, reflection of wake waves is significant. Figure 14(a) shows a time before the ship reaches the depth change centred at $x_a$, so the upstream waves are negligible. As the ship passes the depth variation, mini-tsunamis are induced at the bow (shown in figure 14b), and then propagate ahead of the ship (shown in figure 14c).

Figure 14. The surface elevation $\eta$ at different times in Test 3: (a$t/\sqrt {h_0/g}=139$; (b$t/\sqrt {h_0/g}=186$; (c$t/\sqrt {h_0/g}=224$. Here $h_0$ is the average water depth, and $x_a$ is the centre of the depth change.

Figure 15 shows the free surface along the centreline of the ship region while the ship is sailing through the depth variation. The centre of the depth variation is at $(x-x_a)/L=0$. The ships modelled by the dynamic and static pressure distributions are shown in figures 15(a) and 15(b), respectively. As the ship travels from the deep water to the shallow water, the bottom of the ship modelled by the dynamic pressure stays flat, and the ship volume at different places is nearly the same as the original volume $V_0$, which is shown in figure 16. However, substantial oscillation is shown at the bottom of the ship modelled by the static pressure. Because of the deformation of the hull shape, the mean level of the ship bottom is much deeper than the original ship draft which is at $\eta /d=-1$, inducing a change in the ship volume, which is denoted by the blue line and squares in figure 16. Due to the deformation of the hull shape during the initial acceleration phase, the ship volume is around $V=1.09V_0$ in the relatively deep region. As the ship travels into the shallow-water region, additional deformation occurs, and the ship volume is increased by approximately 27 % compared with that in the relatively deep water. After the ship has passed the depth change and travels in shallow water, the three-dimensional free surface of the ship region modelled by the two types of pressure are compared in figure 17. From these comparisons, the dynamic correction strategy is able to preserve the described hull shape and ensures that the ship has a flat bottom and a stable ship volume.

Figure 15. Surface elevation $\eta$ along the centreline of the ship region as the ship travels through the depth variation. The centre of the depth change is at $(x-x_a)/L=0$. Here $d$ denotes the ship draft. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively.

Figure 16. Variation of the ship volume $V$ as the ship travels over the depth variation. The centre of the depth change is at $(x-x_a)/L=0$. Here $V_0$ is the original ship volume.

Figure 17. Three-dimensional surface $\eta$ in the ship region after the ship has passed the depth change. Here $x_a$ represent the centre of the depth change. (a) Dynamic pressure distribution and (b) static pressure distribution.

During the generation process of the mini-tsunamis, the free surface elevation along the centreline of the channel is shown in figure 18. The wave amplitude of the mini-tsunamis induced by the static pressure distribution is approximately $0.016h_0$. Because the static pressure distribution has an uneven bottom and a larger ship volume than that of the dynamic pressure distribution, the wave amplitude of the mini-tsunamis induced by the dynamic pressure distribution is smaller than that induced by the static pressure distribution. To be specific, the wave amplitude computed by the static pressure distribution is 45 % larger than that obtained using the dynamic pressure distribution. This illustrates that the dynamic correction has a significant impact on the wave amplitude of mini-tsunamis.

Figure 18. Surface elevation $\eta$ along the centreline of the channel as the mini-tsunamis are generated from the bow: (a$t/\sqrt {h_0/g}=182$; (b$t/\sqrt {h_0/g}=190$; (c$t/\sqrt {h_0/g}=207$.

6.2. Mini-tsunami wave characteristics

To systematically study the impact of the dynamic correction on the wave amplitude $A$ and wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis, the test cases from table 1 are considered.

Figure 19 compares the generated wave amplitudes from the two methods as a function of increasing Froude number and bottom slope parameter. At $Fr_0=0.40$, the difference in wave amplitudes induced by the two methods is small, and both results are around $A/h_0=0.005$. As the Froude number increases, the wave amplitude increases, as does the induced fluid acceleration around the ship, and the dynamic effects become increasingly significant. For example, the dynamic correction strategy decreases the wave amplitude by approximately 42 % at $Fr_0=0.55$. This indicates the importance of a dynamic correction when the ship speed is large. Moreover, based on the results for the dynamic pressure and static pressure, fitted curves corresponding to $Fr_0^\alpha$ are obtained, and are denoted by the blue and red lines in figure 19(a), respectively. The value of $\alpha$ is 3.9 for the dynamic pressure distribution and 5.0 for the static pressure distribution. As shown in figure 19(b), the wave amplitude also increases with the slope of the depth variation $\gamma$. However, the decrease in wave amplitude due to the dynamic correction is nearly constant with a value of approximately $0.005h_0$ as the bottom slope changes from 0.7 to 10.15.

Figure 19. Influence of dynamic correction on the amplitudes $A$ of mini-tsunamis in cases with different (a) Froude numbers based on the average water depth $Fr_0$ and (b) slopes of the depth variation $\gamma$.

Figure 20 shows the wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis generated by the dynamic and static pressure distributions. The wavelength diminishes with increasing Froude number and remains nearly constant with varying bottom slope. The dynamic correction results in a slight reduction in wavelength, but the effect is not sensitive to either Froude number or bottom slope. The decrease is less than 6.9 % in all the tested cases.

Figure 20. Influence of dynamic correction on the wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis in cases with different (a) Froude numbers based on the average water depth $Fr_0$ and (b) slopes of the depth variation $\gamma$.

7. Variation of wave resistance when a ship travels through a depth change

The total wave resistance of the ship modelled by the dynamic pressure distribution is computed by the pressure integration method in Tests 1 to 4 with varying Froude numbers based on the average water depth $Fr_0$. The slope of the depth variation $\gamma$ is set to 0.7. The depth variation area along the $x$ direction is shown in figure 21(a). The region delineated by the two blue dashed lines represents a transition zone for the water depth, specifically $50 \leq (x-x_s)/h_0 \leq 70$, where $x_s$ denotes the initial position of the ship. The region defined by $0 \leq (x-x_s)/h_0 < 50$ corresponds to the relatively deep region where the water depth satisfies $h_1/h_0=1.45$, while the region with $50 < (x-x_s)/h_0 < 170$ represents the relatively shallow region where the water depth satisfies $h_2/h_0=0.55$.

Figure 21. (a) Area for the water depth variation along $x$-axis. The red dash–dotted line denotes the centre position of the depth change $x_a$. The transition zone for the water depth is denoted between the two blue dashed lines. The water depths in the deep and shallow regions are $h_1/h_0=1.45$ and $h_2/h_0=0.55$, respectively. $h_0$ is the average water depth. Here $x_s$ is the initial position of the ship. (b) Wave resistance coefficient $C_w$ at different positions and Froude numbers.

The variation in the wave resistance coefficient $C_w$ after the ship initiates motion is depicted in figure 21(b). The ship accelerates from rest and then travels at a constant speed; hence, the wave resistance in the deep region increases from zero to an almost constant value. When the ship traverses the depth change, the wave resistance increases due to the generation of the mini-tsunamis from the bow. However, as the mini-tsunamis detach from the bow and propagate ahead of the ship, the wave resistance decreases, eventually stabilizing. This stable condition is characterized by periodic oscillations of the wave resistance. This oscillation is induced by the short waves behind the mini-tsunamis when the ship travels in the relatively shallow water. Due to the generation and propagation of the mini-tsunamis, a peak of the wave resistance is observed after the ship passes the centre of the depth variation $x_a$. As the Froude number increases, both the peak of the wave resistance and the distance between the position of the resistance peak and the centre of the depth variation also increase. In figure 22, the peak value of $C_w$ is computed using different nonlinear orders $M$. Here $M=1$ means that the linear free-surface boundary condition is used. Second- and third-order nonlinearities increase the peak of the wave resistance, with second-order nonlinearities generally giving the largest contribution. Convergence is shown from results with $M=3$, 4, 5. The nonlinear effects on the wave resistance increase as the Froude number increases. At $Fr_0=0.55$, the resistance peak is increased by more than 56 % due to second- and third-order nonlinearities. This implies that the neglect of the nonlinearity will induce a significant underestimation of wave resistance, especially at large Froude number.

Figure 22. Influence of free-surface nonlinearity on the peak value of the wave resistance coefficient $C_w$ versus Froude number based on the average water depth $Fr_0$. Here $M$ is the nonlinear order of the free surface in the HOS model.

8. Discussion and conclusions

A travelling ship is often modelled by a moving static pressure distribution which neither preserves the desired shape of the hull nor the original volume when the ship moves or encounters a depth change. The ship movement causes changes in the movement of the surrounding flow field, generating dynamic pressure which cannot be balanced by the static pressure approach. By imposing a dynamic correction strategy, we have demonstrated that the hull shape can be accurately maintained under any sailing condition. The influence of the dynamic correction on the ship-induced waves over a flat bottom, upstream mini-tsunamis induced by a ship travelling over an abrupt depth change, and wave resistance is quantified by comparing HOS calculations using the two strategies.

The static pressure approach is theoretically applicable to ACVs, because their hull shape changes according to changes in pressure along the hull surface. The ACV used in the experiments of Doctors (Reference Doctors1993) travels with different Froude numbers, and it is modelled by a moving static pressure distribution. The wave resistance calculated by pressure integration and wave-cut analysis in the present study agrees very well with the experimental data. However, the static pressure approach cannot accurately model ships with rigid hulls.

The in-house HOS code with the dynamic correction strategy is validated by comparing calculations with the experiments of Ertekin et al. (Reference Ertekin, Webster and Wehausen1984), where a Series 60 ship travels along the centreline of a channel with different Froude numbers, channel widths and water depths. In Ertekin et al. (Reference Ertekin, Webster and Wehausen1986), some experimental cases were simulated based on Green–Naghdi equations where a static pressure distribution was proposed to model the ship. In their numerical results, the hull shape is deformed and significant oscillations are shown in the ship region after the ship starts to move, so no precise comparisons were made with experimental values. Here, the dynamic correction strategy, which can accurately calculate the variation of hull-surface pressure induced by the ship movement, is applied to the static pressure distribution. The time series of ship-induced waves and wave resistance obtained using the two types of moving pressure distribution are compared with the experimental records in detail. The upstream waves obtained using the static pressure distribution, dynamic pressure distribution and those in the experiments are comparable. The length of a long depression induced by the hull in the experiments agrees well with that computed by the dynamic pressure distribution, but it is shorter than that computed by the static pressure distribution. The phase of the wake waves using the dynamic strategy compares well with the experiments while the static pressure results are much different. In addition, the amplitude of the wake waves is decreased by the dynamic correction strategy. A correlation coefficient is adopted to quantify the differences between the numerical and experimental results. It is found that the correlation coefficient for the wake waves is improved by approximately 0.25 due to the dynamic correction strategy. Correspondingly, the total wave resistance is also decreased by the dynamic correction strategy, and the resistance induced by the dynamic pressure is more consistent with the experiments, which is approximately 10 % smaller than the resistance induced by the static pressure.

When the ship modelled by a static pressure travels through a depth change, significant oscillation appears along the ship boundary and the ship volume is substantially increased. However, the desired hull shape is well preserved by the dynamic correction strategy. Because of the difference in the ship volume, the dynamic pressure distribution induces mini-tsunamis with a smaller amplitude than those generated by the static pressure distribution. As the Froude number increases, the fluid accelerations induced by the depth-change also increase, and the dynamic correction has an increasingly large impact which can decrease the wave amplitude by around 42 % at $Fr_0=0.55$. Neglecting these changes leads to an overestimation of the amplitudes of ship-induced waves. However, the decrease in the wavelength due to the dynamic correction is not significant, and is smaller than 6.9 % in all the test cases considered here. It should be noted that in our test cases with a variable bottom, the ship travels along the centreline of a narrow channel where the ratio of the channel width to ship width is $W/B=4$. However, Grue (Reference Grue2017) used a much wider channel with $W/B=20.43$, and the static pressure distribution $p(x,y)=\rho{g}d \exp{(-(2x/L)^6-(2y/B)^6)}$ was adopted to simulate the moving ship. The variable topographies in his and the present study are almost the same, but the desired hull shapes are a little different, i.e. (2.9) with dynamic pressure correction is adopted to describe a hull shape in the present study. To investigate the effects of dynamic correction on his results, the same parameters of the numerical simulation are chosen and the comparison is shown in figure 23. The results denoted by red triangles are extracted from figure 9 in Grue (Reference Grue2017). The wave amplitude of mini-tsunamis is slightly reduced by the dynamic correction, with this effect becoming more pronounced as the Froude number increases. These findings are similar to those shown in figure 19(a) with a narrow channel ($W/B=4$), although the differences caused by the dynamic correction are more significant in a narrow channel. In figure 23(b), the reduction in the wavelength due to the dynamic correction is not significant, similar to what is shown in figure 20(a). This indicates that the influence of dynamic correction on mini-tsunamis becomes less significant as the channel width increases.

Figure 23. Influence of dynamic correction on the (a) amplitudes $A$ and (b) wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis in cases studied in Grue (Reference Grue2017), where static pressure distribution was used. The ratio of the channel width to ship width is $W/B=20.43$ and different Froude numbers $Fr_0$ are used.

The wave resistance increases during the generation of a mini-tsunami, then decreases to a new nearly steady condition after the mini-tsunami has separated from the bow. The maximum wave resistance during the generation process is found to be increasingly influenced by free-surface nonlinearity as the Froude number increases. At $Fr_0=0.55$, nonlinear effects increase the peak resistance by more than 56 %. Yao et al. (Reference Yao, Bingham and Zhang2023) emphasized that the amplitude of mini-tsunamis can be substantially increased by nonlinear effects, but they modelled the ship by a static pressure distribution. This paper demonstrates the importance of free-surface nonlinearity on wave resistance when the ship is modelled by a dynamic pressure distribution. Disregarding the impact of free-surface nonlinearity underestimates the wave resistance.

To maintain the desired hull shape for a moving ship, a dynamic correction strategy is introduced by the present study. However, the effects of the sinkage and trim of the ship are not considered here. Due to sinkage and trim, the submerged hull shape differs from that of a stationary ship. In future work, the impact of sinkage and trim will be incorporated into the dynamic correction strategy to enhance the accuracy of simulating ship waves and resistance.

In addition, flow separation at the stern is commonly observed in ships with a transom stern. To ensure the continuity of the free-surface elevation and the pressure along the detaching streamline, Kring (Reference Kring1994) modified the free-surface boundary conditions, and proposed a Kutta condition that include both a kinematic condition and a dynamic condition. The dynamic condition imposes a dynamic pressure on the free surface that balances the hydrostatic pressure resulting from the given transom draft. The total pressure at the separation point is then equal to the atmospheric pressure. Using the Kutta condition, the dynamic correction strategy could be applied to ensure that the pressure along the detaching streamline equals the atmospheric pressure and that the free-surface elevation matches the local draft at the separation point. Subsequently, the effects of flow separation on the free-surface elevation could be investigated using the HOS or other potential-flow models (Kring Reference Kring1994; Mola, Heltai & DeSimone Reference Mola, Heltai and DeSimone2017).

Acknowledgements

The authors appreciate the valuable discussions with Professor R.C. Ertekin at the University of Hawaii.

Funding

X.Z. and J.Y. gratefully acknowledge the financial support from the National Natural Science Foundation of China under grant no. 52171269.

Declaration of interests

The authors report no conflict of interest.

Appendix A. Savitzky–Golay filter

Suppose there is a string of data with a frame length of $L_f=2k_f+1$ to be smoothed, namely, $\boldsymbol {Y}=[y_{j-k_f},\ldots,y_{j-1},y_j,y_{j+1},\ldots,y_{j+k_f}]^{\rm T}$, where $k_f$ is a positive integer and the superscript T means the transpose of a matrix or vector. These data are uniformly distributed along the $x$-axis. The corresponding positions along the $x$-axis are $\boldsymbol {X}=[x_{j-k_f},\ldots,x_{j-1},x_j,x_{j+1},\ldots,x_{j+k_f}]^{\rm T}$, and the distance between two adjacent points is $\Delta {x}$. Based on the method of a SG filter, the estimation of this string of data is calculated by least-squares fitting an $n_f$th-order polynomial through the values of data in the frame length. Note that the frame length $L_f$ is a positive odd integer greater than the order $n_f$. Therefore, this string of data is expressed by

(A1)\begin{equation} \left[ \begin{array}{@{}c@{}} y_{j-k_f} \\ \vdots \\ y_{j-1} \\ y_j \\ y_{j+1} \\ \vdots \\ y_{j+k_f} \\ \end{array} \right] = \left[\begin{array}{@{}c@{}} \displaystyle\sum_{i=0}^{n_f} b_i(x_j-k_f\Delta{x})^i \\ \vdots \\ \displaystyle\sum_{i=0}^{n_f} b_i(x_j-1\Delta{x})^i \\ \displaystyle\sum_{i=0}^{n_f} b_i(x_j-0\Delta{x})^i \\ \displaystyle\sum_{i=0}^{n_f} b_i(x_j+1\Delta{x})^i \\ \vdots \\ \displaystyle\sum_{i=0}^{n_f} b_i(x_j+k_f\Delta{x})^i \\ \end{array} \right] = \left[\begin{array}{@{}c@{}} \displaystyle\sum_{i=0}^{n_f} a_i(-k_f)^i \\ \vdots \\ \displaystyle\sum_{i=0}^{n_f} a_i(-1)^i \\ \displaystyle\sum_{i=0}^{n_f} a_i(0)^i \\ \displaystyle\sum_{i=0}^{n_f} a_i(1)^i \\ \vdots \\ \displaystyle\sum_{i=0}^{n_f} a_i(k_f)^i \end{array} \right], \end{equation}

where $b_i$ is the unknown coefficient corresponding to each term $(x_j+l_f\Delta {x})^i$ in the polynomial $\sum _{i=0}^{n_f} b_i(x_j+l_f\Delta {x})^i$, and $l_f=-k_f, -(k_f-1), \ldots, -1, 0, 1, \ldots, k_f-1, k_f$. This polynomial can be expanded as $\sum _{i=0}^{n_f} [b_i\sum _{r=0}^{i}{\mathcal {C}}_{i}^r(x_j)^{i-r}(l_f\Delta {x})^r]$, where $\mathcal {C}$ is the combination formula. Changing the basis $(x_j-l_f\Delta {x})^i$ to $(l_f)^i$, the polynomial $\sum _{i=0}^{n_f} b_i(x_j+l_f\Delta {x})^i$ is rewritten as $\sum _{i=0}^{n_f} a_i(l_f)^i$. The coefficient $a_i=\sum _{s=i}^{n_f} b_{s}{\mathcal {C}}_{s}^i(x_i)^{s-i}(\Delta {x})^i$. Based on (A1), the original data string $\boldsymbol {Y}$ is written as

(A2)\begin{equation} \boldsymbol{Y} = \left[\begin{array}{@{}ccccc@{}} 1 & -k_f & (-k_f)^2 & \cdots & (-k_f)^{n_f} \\ 1 & \cdots & \cdots & \ddots & \vdots \\ 1 & -2 & (-2)^2 & \cdots & (-2)^{n_f}\\ 1 & -1 & (-1)^2 & \cdots & (-1)^{n_f} \\ 1 & 0 & 0 & \cdots & 0 \\ 1 & 1 & (1)^2 & \cdots & (1)^{n_f} \\ 1 & 2 & 2^2 & \cdots & 2^{n_f}\\ 1 & \cdots & \cdots & \ddots & \vdots \\ 1 &k_f & k_f^2 & \cdots & k_f^{n_f} \\ \end{array} \right] \left[ \begin{array}{@{}c@{}} a_0 \\ \vdots \\ a_n \end{array} \right] = \boldsymbol{HA}, \end{equation}

where the coefficient vector $\boldsymbol {A}=[a_0, a_1, \ldots, a_n]^{\rm T}$. To find the SG estimates of the original data, use the pseudoinverse of $\boldsymbol {H}$ to compute the coefficient vector $\boldsymbol {A}$ and then premultiply by $\boldsymbol {H}$. Therefore, the estimation of the original data is solved by

(A3)\begin{equation} \boldsymbol{\hat{Y}}=\boldsymbol{H} \left( \boldsymbol{H}^{\rm T}\boldsymbol{H} \right)^{-1}\boldsymbol{H}^{\rm T}\boldsymbol{Y}, \end{equation}

where $\boldsymbol {\hat {Y}}$ is the results obtained by applying a SG filter on the original data $\boldsymbol {Y}$.

To mitigate the sawtooth instabilities in the HOS model, a SG filter with an order of $n_f$ as large as possible is used. Increasing the order of $n_f$ results in less energy being filtered out. As shown in figure 24, the order $n_f$ gradually increases with the distance from the ship region, and the frame length is a constant, i.e. $L_f=11$. Specifically, at each time step in the HOS model with a dynamic or static pressure distribution, each row of nodes along the $x$-axis on the free surface is first filtered by the SG filter shown in figure 24(a), and then each column of nodes along the $y$-axis is filtered by the SG filter shown in figure 24(b).

Figure 24. The order $n_f$ of the SG filter in the whole computational domain. $L$ and $B$ are the ship length and width, respectively. The red square represents the centre of the water plane of the ship at time $t$. (a) The orders $n_f$ of the SG filters which are applied on each row of nodes along the $x$-axis in the computational domain. (b) The orders $n_f$ of the SG filters which are applied on each column of the nodes along the $y$-axis in the computational domain.

Appendix B. The influence of the penalty term on the hull shape of a moving ship

The last term in (2.14) is a penalty term, i.e. $\rho {g}(\eta -\eta _s)$, to help enforce the condition that $\eta =\eta _s$ in the ship region and improve the stability of the dynamic correction strategy. We recommend including this penalty term, especially when the Froude number is large. To illustrate the effect of this term, we multiply it by a weighting factor $\lambda _{pen}$, i.e. $\lambda _{pen}\rho {g}(\eta -\eta _s)$. If $\lambda _{pen}=0$, no penalty term is added. If $\lambda _{pen}=1$, the penalty term is at fully strength. The influence of the penalty term on the hull shape is shown in figure 25 where the ship adopted in the experiments of Ertekin et al. (Reference Ertekin, Webster and Wehausen1984) is used as an example. The detailed information about the ship and experiments can be found in § 3.1. The ship travels along the centreline of the narrow channel and the Froude number $Fr=1.1$. The magenta line exhibits the desired hull shape. Different weighting factors, i.e. $\lambda _{pen}=0$, 0.5 and 1, are used in the numerical simulations, respectively. The hull shape along the centreline of the ship region at $t=21\sqrt {(h/g)}$ is shown in figure 25. With the penalty term at full strength, the desired hull shape is well preserved while the ship is moving. However, if the penalty term is discarded, i.e. $\lambda _{pen}=0$, the ship bottom shows some slight oscillations which has an impact on the stability of the dynamic correction strategy, although the averaged ship draft is close to the initial draft. If half of the penalty term is added, the slight oscillation at the ship bottom is reduced, but the result is still not as good as with the term at full strength. Therefore, it is recommended to include the penalty term in the numerical simulation.

Figure 25. The surface elevation $\eta$ along the centreline of the ship region. The magenta line represents the desired hull shape at the initial time when the ship is stationary. Different values of a weighting factor $\lambda _{pen}$ are used in the dynamic correction strategy to simulate a moving ship. The red, black and blue lines illustrate the travelling ship at $t=21\sqrt {(h/g)}$ with varying weighting factors. Here $d$, $x_s$ and $h$ denote the ship's draft, the initial position of the ship and the water depth over a flat topography, respectively.

Appendix C. Methods for wave resistance calculation

C.1. Direct pressure integration

The total wave resistance can be directly calculated by a pressure integration along wetted hull surface,

(C1)\begin{equation} R_w=\iint_{S_B} p(x,y,z)n_x \,{\rm d}S, \end{equation}

where $S_B$ is the wetted hull surface, $p$ is the pressure along the wetted hull surface and $n_x$ is the component of the normal vector along the $x$-axis. In the present study, $S_B$ is strictly confined to the ship region where the ship pressure is applied. Consequently, the resistance induced by the wake waves and the upstream waves cannot be distinguished by the pressure integration method.

C.2. Far-field method

Consider a control volume $V^c$ bounded by the wetted hull surface $S_B$, the free surface $S_F$ and the far-field control surface $S_{\infty }$. To be specific, $S_{\infty }$ includes the sea bottom plane $S_{S}$, a plane ahead of the bow $S_A$, a plane behind the stern $S_{W}$ and two sidewalls of the tank $S_T$. The far-field method follows the momentum conservation principle (Kashiwagi Reference Kashiwagi2022). In one period, the time-averaged value of change of momentum $Q$ in the control volume $V^c$ is zero,

(C2)\begin{equation} \overline{\frac{{\rm d}Q_i}{{\rm d}t}}=\overline{\iint_{S_B+S_F+S_{\infty}} \left\{ pn_i+\rho{v}_i({v}_n-u_n) \right\} \,{\rm d}S }=0, \end{equation}

where the overbar denotes the time average, $Q_i$ is the momentum in the $i$-direction, which can correspond to the $x$-axis or the $y$-axis. Similarly, $n_{i}$ represents the component of the normal vector in the $i$-direction, $v$ represents the fluid velocity, $u_n$ represents the normal velocity of the boundary surface. Therefore, the total wave resistance is expressed by

(C3)\begin{align} \overline{F_x}&=\overline{ \iint_{S_B} p{n}_x \,{\rm d}S }\nonumber\\ &=-\overline{ \iint_{S_{\infty}} \left\{ p{n}_x+\rho\phi_x(\phi_n-u{n}_x) \right\}\, {\rm d}S }.\end{align}

The time-averaged value of this integration on the plane ahead of the bow $S_A$ is used as the resistance induced by the upstream waves.

C.3. Wave-cut analysis

In the reference frame fixed to the ship, wake waves generated by a ship moving with a steady speed $U$ are assumed to be stationary (Wehausen Reference Wehausen1973). Within a channel of width $W$, the free surface elevation $\eta$ and the velocity potential $\phi$ are expressed as follows:

(C4)$$\begin{gather} \eta(x,y)=\sum_{n=0}^{\infty}\left[a_n\sin(k_nx)+b_n\cos(k_nx)\right]\cos\frac{n{\rm \pi}}{W}\left(y-\frac{W}{2}\right), \end{gather}$$
(C5)$$\begin{gather}\phi(x,y)=\sum_{n=0}^{\infty}\frac{g}{Uk_n}\frac{1}{\cosh{\mu_nh}}\nonumber\\ \times\left[-a_n\cos(k_nx)+b_n\sin(k_nx)\right]\cosh\left[\mu_n(z+h)\right]\cos\frac{n{\rm \pi}}{W}\left(y-\frac{W}{2}\right), \end{gather}$$

where the two coefficients $a_n$ and $b_n$ are related to the wave amplitude, $h$ is the water depth. The wavenumber $k_n$ satisfies $k_n^2=\mu _n^2-(n{\rm \pi} /W)^2$, and $\mu _n$ is calculated by

(C6)\begin{equation} \frac{U^2}{gh}\left[\mu_nh-\left(\frac{n{\rm \pi}{h}}{W}\right)^2\frac{1}{\mu_nh}\right]=\tanh{\mu_nh}. \end{equation}

The wave-cut analysis is based on the momentum conservation principle. Assuming that wake waves are stationary in the ship-fixed reference frame and no upstream waves are generated, the expression of the resistance in (C3) simplifies to

(C7)\begin{equation} F_x=\tfrac{1}{2}\rho \iint_{S_{W}} \left(-\phi_x^2+\phi_y^2+\phi_z^2\right) \,{\rm d}S+ \tfrac{1}{2}\rho{g} \int_{-W/2}^{W/2} \eta(x_W,y) \,{\rm d} y, \end{equation}

where $x_W$ is the position of the wave-cut plane $S_W$ along the $x$-direction. By substituting (C4) and (C5) into (C7), the resistance induced by wake waves is

(C8)\begin{equation} R_{ww}=\frac{1}{4}\rho{g}W\sum_{n=0}^{\infty}\varepsilon_n^{-1}\left(a_n^2+b_n^2\right)\left[1-\frac{1}{2}\frac{gh}{U^2}\frac{\tanh\mu_nh}{\mu_nh}\left(1+\frac{2\mu_nh}{\sinh2\mu_nh}\right) \right], \end{equation}

where $\varepsilon _0=\frac {1}{2}$, $\varepsilon _n=1$ for $n\geq 1$.

Appendix D. Illustration about equivalent channel walls

In figure 26, the area enclosed by the black solid lines represents the computational domain used in the numerical simulation. The length and width of the computational domain are $L_x$ and $L_y$, respectively. The width $L_y$ equals the channel width $W$ in each test case. Therefore, the computational domain includes one channel, through which a ship travels along the centreline. Since periodic boundary conditions are used in the horizontal plane in the HOS model, the computational domain can be periodically shown along the $y$-axis. Consequently, virtual domains bounded by the black dashed lines are constructed based on this periodic boundary condition. The size of each virtual domain is the same as the computational domain. The two virtual domains situated on the left- and right-hand sides of the ship within the computational domain are denoted as $1_L$ and $1_R$, respectively. Theoretically, an infinite number of such virtual domains can exist, with $1_L$ and $1_R$ being the closest two domains to the computational domain. The virtual ship, indicated by the blue dashed lines, travels along the centreline within each virtual domain. The trajectories of two adjacent ships are symmetric with respect to the boundary positioned between them in the middle, meaning the paths of the virtual ship in virtual domain $1_L$ or $1_R$ mirror those of the ship in the computational domain across the boundary along the $x$-axis between them. Hence, the two boundaries of the computational domain along the $x$-axis can be regarded as the channel walls.

Figure 26. Illustration about equivalent channel walls. The area enclosed by the black solid lines represents the computational domain, with dimensions $L_x$ in length and $L_y$ in width. Virtual domains, bounded by black dashed lines, are created according to the periodic boundary conditions. Within each virtual domain, the virtual ship, marked by blue dashed lines, travels along the centreline. The trajectories of the ships are depicted by red dot–dashed lines.

References

Amini-Afshar, M. & Bingham, H.B. 2021 Added resistance using Salvesen–Tuck–Faltinsen strip theory and the Kochin function. Ocean Engng 106, 102481.Google Scholar
Amini-Afshar, M., Bingham, H.B. & Kashiwagi, M. 2023 On a new formulation for wave added resistance. In Proceedings of 38th International Workshop on Water Waves and Floating Bodies, Ann Arbor, USA.Google Scholar
Dias, F. 2014 Ship waves and Kelvin. J. Fluid Mech. 746, 14.CrossRefGoogle Scholar
Doctors, L.J. 1972 The wave resistance of an air-cushion vehicle in steady and accelerated motion. J. Ship Res. 16, 248260.CrossRefGoogle Scholar
Doctors, L.J. 1993 On the use of pressure distributions to model the hydrodynamics of air-cushion vehicles and surface-effect ships. Nav. Engng J. 105 (2), 6989.CrossRefGoogle Scholar
Dommermuth, D.G. & Yue, D.K.P. 1987 A high-order spectral method for the study of nonlinear gravity waves. J. Fluid Mech. 184, 267288.CrossRefGoogle Scholar
Ertekin, R.C., Webster, W.C. & Wehausen, J.V. 1984 Ship-generated solitons. In Proceedings of the 15th Symposium on Naval Hydrodynamics, Hamburg, Germany, pp. 22–36. National Academy of Sciences.Google Scholar
Ertekin, R.C., Webster, W.C. & Wehausen, J.V. 1986 Waves caused by a moving disturbance in a shallow channel of finite width. J. Fluid Mech. 169, 275292.CrossRefGoogle Scholar
Gouin, M., Ducrozet, G. & Ferrant, P. 2017 Propagation of 3D nonlinear waves over an elliptical mound with a high-order spectral method. Eur. J. Mech. (B/Fluids) 63, 924.CrossRefGoogle Scholar
Gourlay, T.P. 2001 The supercritical bore produced by a high-speed ship in a channel. J. Fluid Mech. 434, 399409.CrossRefGoogle Scholar
Grue, J. 2017 Ship generated mini-tsunamis. J. Fluid Mech. 816, 142166.CrossRefGoogle Scholar
Grue, J. 2020 Mini-Tsunami made by ship moving across a depth change. ASCE J. Waterway Port Coastal Ocean Engng 146 (5), 04020023.CrossRefGoogle Scholar
Grue, J., Pedersen, G.K. & Saetra, Ø. 2022 Free wave effects in meteotsunamis. J. Geophys. Res. 127, e2021JC017669.CrossRefGoogle Scholar
Huang, T.T. & Wong, K.K. 1970 Disturbance induced by a pressure distribution moving over a free surface. J. Ship Res. 14 (3), 195203.CrossRefGoogle Scholar
Janson, C.E. & Spinney, D.A. 2004 Comparison of four wave cut analysis methods for wave resistance prediction. Ship Technol. Res. 51, 173184.CrossRefGoogle Scholar
Kashiwagi, M. 2022 A new computation method for added resistance and connection with Lagally's theorem. In Proceedings of the Japan Society of Naval Architects and Ocean Engineers, Kobe, Japan, pp. 253–259. The Japan Society of Naval Architects and Ocean Engineers.Google Scholar
Kelvin, Lord 1887 On ship waves. In Proceedings of the Institution of Mechanical Engineers, Edinburgh, UK, pp. 409–434. Institution of Mechanical Engineers.CrossRefGoogle Scholar
Kring, D.C. 1994 Time domain ship motions by a three-dimensional Rankine panel method. PhD thesis, Massachusetts Institute of Technology, Cambridge, MA, USA.Google Scholar
Li, Y. & Sclavounos, P.D. 2002 Three-dimensional nonlinear solitary waves in shallow water generated by an advancing disturbance. J. Fluid Mech. 470, 383410.CrossRefGoogle Scholar
Lindberg, O., Glimberg, S.L., Bingham, H.B., Engsig-Karup, A.P. & Schjeldahl, P.J. 2013 Towards real time simulation of ship-ship interaction - Part II: Double body flow linearization and GPU implementation. In Proceedings of the 28th International Workshop on Water Waves and Floating Bodies, L'Isle sur la Sorgue, France.Google Scholar
Liu, S., Waseda, T., Yao, J. & Zhang, X. 2022 a Statistical properties of surface gravity waves and freak wave occurrence in crossing sea states. Phys. Rev. Fluids 7, 074805.CrossRefGoogle Scholar
Liu, S., Waseda, T. & Zhang, X. 2022 b Four-wave resonant interaction of surface gravity waves in finite water depth. Phys. Rev. Fluids 7, 114803.CrossRefGoogle Scholar
Ma, C., Zhang, C., Huang, F., Yang, C., Gu, X., Li, W. & Noblesse, F. 2017 Practical evaluation of sinkage and trim effects on the drag of a common generic freely floating monohull ship. Appl. Ocean Res. 65, 111.CrossRefGoogle Scholar
Mola, A., Heltai, L. & DeSimone, A. 2017 Wet and dry transom stern treatment for unsteady and nonlinear potential flow model for naval hydrodynamics simulations. J. Ship Res. 61 (1), 114.CrossRefGoogle Scholar
Noblesse, F., He, J., Zhu, Y., Liang, H., Zhang, C., Zhu, R. & Yang, C. 2014 Why can ship wakes appear narrower than Kelvin's angle? Eur. J. Mech. (B/Fluids) 46, 164171.CrossRefGoogle Scholar
Pethiyagoda, R., McCue, S.W. & Moroney, T.J. 2015 Wake angle for surface gravity waves on a finite depth fluid. Phys. Fluids 27, 061701.CrossRefGoogle Scholar
Pethiyagoda, R., Moroney, T.J., Lustri, C.J. & McCue, S.W. 2021 Kelvin wake pattern at small Froude numbers. J. Fluid Mech. 915, A126.CrossRefGoogle Scholar
Shi, F., Malej, M., Smith, J.M. & Kirby, J.T. 2018 Breaking of ship bores in a Boussinesq-type ship-wake model. Coastal Engng 132, 112.CrossRefGoogle Scholar
Tanimoto, K., Kobayashi, H. & Ca, V.T. 2000 Ship waves in a shallow and narrow channel. In Proceedings of the 27th International Conference on Coastal Engineering, Sydney, Australia (ed. B.L. Edge), pp. 1141–1154. American Society of Civil Engineers.CrossRefGoogle Scholar
Wehausen, J.V. 1973 The wave resistance of ships. Adv. Appl. Mech. 13, 93245.CrossRefGoogle Scholar
West, B., Brueckner, K., Janda, R.S., Milder, D.M. & Milton, R.L. 1987 A new numerical method for the surface hydrodynamics. J. Geophys. Res. 92, 1180311824.CrossRefGoogle Scholar
Wu, T.Y. 1987 Generation of upstream advancing solitons by moving disturbances. J. Fluid Mech. 184, 7599.CrossRefGoogle Scholar
Wu, H., He, J., Liang, H. & Noblesse, F. 2019 Influence of Froude number and submergence depth on wave patterns. Eur. J. Mech. (B/Fluids) 75, 258270.CrossRefGoogle Scholar
Yao, J., Bingham, H.B. & Zhang, X. 2023 Nonlinear effects of variable bathymetry and free surface on mini-tsunamis generated by a moving ship. Phys. Rev. Fluids 8, 094801.CrossRefGoogle Scholar
Yeung, R.W., Wan, H., Banumurthy, S.P., Ham, W.L. & Lew, J.M. 2008 Effects of configuration on wave resistance of multiple pressure distributions. Mar. Syst. Ocean Technol. 4 (2), 5362.CrossRefGoogle Scholar
Zhu, Y., He, J., Zhang, C., Wu, H., Wan, D., Zhu, R. & Noblesse, F. 2015 Farfield waves created by a monohull ship in shallow water. Eur. J. Mech. (B/Fluids) 49, 226234.CrossRefGoogle Scholar
Figure 0

Figure 1. Sketch of the experiments on solitons over a flat topography in the (a) narrow and (b) wide channels (Ertekin et al.1984). The ship moves along the centreline of the channel with a speed $U{(t)}$. The blue lines indicate the ship. Here $L$ and $B$ denote the length and width of the ship, respectively. Red circles with crosses at the centre represent the WGs.

Figure 1

Figure 2. Sketch of numerical simulation for the experiments conducted by Ertekin et al. (1984). The ship moves along the centreline of the channel at a speed $U{(t)}$. The blue lines indicate the ship. Here $L$ and $B$ denote the length and width of the ship, respectively. The four red circles with crosses at the centre indicate the positions of the four NWGs.

Figure 2

Figure 3. The positions of the wave-cut plane and the far-field plane used for calculation of wave resistance. Here $L$ is the ship length, and $U$ is the ship velocity.

Figure 3

Figure 4. Sketch of the test cases on mini-tsunamis over a variable topography. (a) Plan view: the ship moves along the centreline of the channel. The red square represents the centre of the water plane, whose coordinates are $(x^*,y^*)$. The ship accelerates in the relatively deep region, and moves at a constant velocity through a transition in water depth towards the shallower region. The regions between the two blue dashed lines are the transition zone for the water depth. Here $x_a$ and $x_b$ are the centre positions of the two depth changes. (b) Side view: $d$ is the ship draft. The initial position of the ship is denoted by $x_s$. Here $h_1$ and $h_2$ are the depths in the deep and shallow regions, respectively.

Figure 4

Table 1. Main parameters used in the variable bathymetry test cases. Here $Fr_0$ denotes the Froude number based on the average depth and $\gamma$ controls the slope of the depth variation.

Figure 5

Figure 5. Wave resistance coefficient $C_w=\rho {g}R_w/p_0^2B$ for an ACV as a function of Froude number $Fr_L=U/\sqrt {gL}$. The ACV has $L/B=3/2$ and travels in open water with $h/L=3.58$. The red circles indicate pressure integration results while the blue crosses show the wave-cut analysis results. The lines are splines through the point values. The black triangles show the experimental data in Doctors (1993).

Figure 6

Figure 6. Comparison of the wave gauge records computed by the HOS ($M=3$) with the dynamic correction strategy and those from experiments. The water depth satisfies $h/d=1.67$. The Froude numbers $Fr=U/\sqrt {gh}$ are (a) $Fr=0.8$, (b) $Fr=1.0$ and (c) $Fr=1.1$. The blue solid lines are computed results and the black circles are the experimental data from Ertekin et al. (1984).

Figure 7

Figure 7. The variation of the non-dimensional amplitude $A/h$ of the main upstream wave with the Froude number. The water depth is $h/d=2$. The black circles represent the experimental records from Ertekin et al. (1984). The blue squares represent the numerical results computed by HOS ($M=3$) with the dynamic correction strategy.

Figure 8

Figure 8. The surface elevation $\eta$ along the centreline of the ship region at the initial stage of ship movement. The initial position of the ship is at $(x-x_s)/L=0$. Here $d$ denotes the ship draft. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively.

Figure 9

Figure 9. Variation of the ship volume $V$ at the initial stage of ship movement. The initial position of the ship is at $(x-x_s)/L=0$. Here $V_0$ is the original ship volume.

Figure 10

Figure 10. Time series of free surface elevation $\eta$ recorded by WG3. The whole gauge records comprise the upstream waves, the depression and the wake waves. The black dash–dotted line represents the experimental data from Ertekin et al. (1984). The blue and red dashed lines denote the results computed by the dynamic pressure and static pressure, respectively. Here (a) $Fr=1.0$, the depression is between $195\leq {t/\sqrt {h/g}}\leq 221$; (b) $Fr=1.1$, the depression is between $186\leq {t/\sqrt {h/g}}\leq 204$.

Figure 11

Table 2. Correlation coefficients $\epsilon$ between the numerical and experimental results shown in figure 10. The correlation coefficients of upstream waves, wake waves and the whole gauge records for the dynamic and static pressure distributions are calculated, respectively.

Figure 12

Figure 11. The free surface elevation $\eta$ along the centreline of the channel with a width of $W/B=5.21$ at $t=133\sqrt {h/g}$. The water depth satisfies $h/d=2$.

Figure 13

Figure 12. The variation of the instantaneous resistance coefficient generated by the upstream solitons $C'_{uw} (t)$ with time. Here $T_s$ is the period of the solitons. The channel width and the water depth satisfy $W/B=5.21$ and $h/d=2$, respectively.

Figure 14

Figure 13. Variation of the wave resistance coefficients over one soliton period $T_s$. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively. The blue solid line represents the resistance coefficient from the wake waves $C_{ww}$. The red solid line represents the coefficient of the instantaneous resistance induced by the upstream waves $C'_{uw} (t)$, and its time-averaged value in one period $C_{uw}=\overline {C'_{uw} (t)}$ is denoted by the red dashed line. The coefficient of the instantaneous total wave resistance computed by $C'_{uw} (t)+C_{ww}$ is denoted by the black solid line, its time-averaged value in one period $C_{w}=\overline {C'_{uw} (t)+C_{ww}}$ is denoted by the black dashed line. The magenta dashed line and the green dash–dotted line represent the total wave resistance coefficients $C_w$ obtained by the pressure integration (PI) method and the experiments (EXP), respectively.

Figure 15

Table 3. Comparison of the resistance coefficients computed by the numerical models and those in the experiments. Here $C_{ww}$ and $C_{uw}$ are the coefficients of resistance induced by the wake waves and upstream waves, respectively. Here $C_w$ is the total wave resistance coefficient. Calculation methods are indicated in parentheses. Here ‘WA’, ‘FM’, ‘PI’ and ‘EXP’ denote ‘wave-cut analysis’, ‘far-field method’, ‘pressure integration’ and ‘experiment’, respectively.

Figure 16

Figure 14. The surface elevation $\eta$ at different times in Test 3: (a$t/\sqrt {h_0/g}=139$; (b$t/\sqrt {h_0/g}=186$; (c$t/\sqrt {h_0/g}=224$. Here $h_0$ is the average water depth, and $x_a$ is the centre of the depth change.

Figure 17

Figure 15. Surface elevation $\eta$ along the centreline of the ship region as the ship travels through the depth variation. The centre of the depth change is at $(x-x_a)/L=0$. Here $d$ denotes the ship draft. The ships are modelled by (a) the dynamic pressure distribution and (b) the static pressure distribution, respectively.

Figure 18

Figure 16. Variation of the ship volume $V$ as the ship travels over the depth variation. The centre of the depth change is at $(x-x_a)/L=0$. Here $V_0$ is the original ship volume.

Figure 19

Figure 17. Three-dimensional surface $\eta$ in the ship region after the ship has passed the depth change. Here $x_a$ represent the centre of the depth change. (a) Dynamic pressure distribution and (b) static pressure distribution.

Figure 20

Figure 18. Surface elevation $\eta$ along the centreline of the channel as the mini-tsunamis are generated from the bow: (a$t/\sqrt {h_0/g}=182$; (b$t/\sqrt {h_0/g}=190$; (c$t/\sqrt {h_0/g}=207$.

Figure 21

Figure 19. Influence of dynamic correction on the amplitudes $A$ of mini-tsunamis in cases with different (a) Froude numbers based on the average water depth $Fr_0$ and (b) slopes of the depth variation $\gamma$.

Figure 22

Figure 20. Influence of dynamic correction on the wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis in cases with different (a) Froude numbers based on the average water depth $Fr_0$ and (b) slopes of the depth variation $\gamma$.

Figure 23

Figure 21. (a) Area for the water depth variation along $x$-axis. The red dash–dotted line denotes the centre position of the depth change $x_a$. The transition zone for the water depth is denoted between the two blue dashed lines. The water depths in the deep and shallow regions are $h_1/h_0=1.45$ and $h_2/h_0=0.55$, respectively. $h_0$ is the average water depth. Here $x_s$ is the initial position of the ship. (b) Wave resistance coefficient $C_w$ at different positions and Froude numbers.

Figure 24

Figure 22. Influence of free-surface nonlinearity on the peak value of the wave resistance coefficient $C_w$ versus Froude number based on the average water depth $Fr_0$. Here $M$ is the nonlinear order of the free surface in the HOS model.

Figure 25

Figure 23. Influence of dynamic correction on the (a) amplitudes $A$ and (b) wavelength at half-crest height $\lambda _{1/2}$ of mini-tsunamis in cases studied in Grue (2017), where static pressure distribution was used. The ratio of the channel width to ship width is $W/B=20.43$ and different Froude numbers $Fr_0$ are used.

Figure 26

Figure 24. The order $n_f$ of the SG filter in the whole computational domain. $L$ and $B$ are the ship length and width, respectively. The red square represents the centre of the water plane of the ship at time $t$. (a) The orders $n_f$ of the SG filters which are applied on each row of nodes along the $x$-axis in the computational domain. (b) The orders $n_f$ of the SG filters which are applied on each column of the nodes along the $y$-axis in the computational domain.

Figure 27

Figure 25. The surface elevation $\eta$ along the centreline of the ship region. The magenta line represents the desired hull shape at the initial time when the ship is stationary. Different values of a weighting factor $\lambda _{pen}$ are used in the dynamic correction strategy to simulate a moving ship. The red, black and blue lines illustrate the travelling ship at $t=21\sqrt {(h/g)}$ with varying weighting factors. Here $d$, $x_s$ and $h$ denote the ship's draft, the initial position of the ship and the water depth over a flat topography, respectively.

Figure 28

Figure 26. Illustration about equivalent channel walls. The area enclosed by the black solid lines represents the computational domain, with dimensions $L_x$ in length and $L_y$ in width. Virtual domains, bounded by black dashed lines, are created according to the periodic boundary conditions. Within each virtual domain, the virtual ship, marked by blue dashed lines, travels along the centreline. The trajectories of the ships are depicted by red dot–dashed lines.