A numerical investigation of momentum flux and kinetic energy transfers between turbulent wind and propagating waves

Abstract

This paper focuses on simulating turbulent flow over propagating waves by solving the full Navier–Stokes equations in a moving frame. A careful comparison of flow statistics with previous experimental and numerical results demonstrates, to some extent, the rationality of simplifying wind waves as turbulent flow over moving wave boundaries. The phase-averaging method is then applied to investigate the momentum and energy transfers between turbulent wind and waves propagating at slow, intermediate and fast speeds. The results suggest that the dominant mechanism for producing Reynolds shear stress (RSS) and turbulent kinetic energy (TKE) is related to the wave age. Slow waves produce RSS and TKE similar to a two-dimensional shear turbulence. However, a fast wave enhances the streamwise Reynolds normal stress, the windward side's negative RSS and the gradient of both streamwise and vertical velocities, leading to additional RSS and TKE productions that can be ignored under the slow wave regimes. A strengthening wave–turbulence exchange is also found for fast waves. The intermediate wave can be regarded as a transitional condition determining this change.

Impact Statement

The wind–wave problem has long been an open question due to the elusiveness of its strongly coupled dynamics, hindering our deeper understanding of environmental complex interface flow. This research advances our understanding of turbulent flow over propagating waves through a comprehensive simulation approach employing the full Navier–Stokes equations in a moving frame. Through validations against both experimental and numerical data, the study rationalizes the representation of wind waves as turbulent flow over moving wave boundaries to a certain extent. The findings elucidate distinct wind–wave dynamical behaviours in Reynolds shear stress production and kinetic energy transfer. The study also reveals the different dominant mechanisms dictating momentum flux and kinetic energy productions, contingent upon a transitional condition: intermediate wave.

1. Introduction

The interaction between turbulent wind and water surface waves plays a significant role in shaping the characteristics of both wind and waves. This dynamic process primarily occurs near the air–sea interface and governs the exchange of mass, momentum and energy between the upper ocean and lower atmosphere. Therefore, gaining a comprehensive understanding of these dynamics, particularly regarding momentum and energy transfers, is essential for enhancing our insight into the intricate interactions between wind and waves.

Wind waves can be classified as slow, intermediate and fast waves based on their wave ages (characterized by the ratio of the surface wave's phase velocity $c$ to the friction velocity $u_{\tau }$), corresponding to ${c}/{{{u}_{\tau }}}\le 15$, ${c}/{{{u}_{\tau }}}\approx 15$ and ${c}/{{{u}_{\tau }}}\ge 15$, respectively. Due to the universality of growing wind waves, the slow wave regime has received significant attention in understanding momentum transfer. Experimental (Snyder et al. 1981; Fairall et al. 1996; Hristov, Miller & Friehe 2003; Donelan et al. 2006; Babanin et al. 2007; Longo & Losada 2012; Buckley & Veron 2016; Vollestad, Ayati & Jensen 2019) and numerical (Sullivan, McWilliams & Moeng 2000; Yang & Shen 2009; Druzhinin, Troitskaya & Zilitinkevich 2012; Akbarzadeh & Borazjani 2019; Åkervik & Vartdal 2019) methods were employed to investigate how momentum is transferred in this regime. In general, a slow wave extracts momentum from turbulent wind through the interaction of wave-correlated pressure against the surface-wave slope (Åkervik & Vartdal 2019). This process is accompanied by variations in turbulent stress and wave-induced stress. Belcher & Hunt (1993) discovered that the turbulent stress exhibits wave-coherent characteristics, with positive–negative variation along the wave. This behaviour is similar to what is observed in stationary wavy wall turbulence (Hudson, Dykhno & Hanratty 1996; Cherukat et al. 1998; Yang & Shen 2010; Hamed et al. 2015; Zhang, Wang & Liu 2022a; Zhang et al. 2022b). Consequently, it leads to the typical asymmetrical flow acceleration over the wave, which creates favourable conditions for wave growth (Sullivan & McWilliams 2010). Additionally, studies were carried out to examine the changes in wave-induced stress at the critical height, where the wave-coherent velocity matches the phase velocity (Sullivan et al. 2000; Hristov et al. 2003; Grare, Lenain & Melville 2013a). It was observed that, in the case of slow waves, the wave-induced stress beneath the critical height is primarily positive, but it decreases and becomes negative above the critical height (Yousefi, Veron & Buckley 2020a). However, fast waves exhibit a different distribution of turbulent and wave-induced stresses. Yang & Shen (2010) numerically showed that fast waves result in a more symmetric distribution of turbulent stress. Furthermore, Hristov et al. (2003) observed a symmetrical wave-induced stress pattern above a fast wave, suggesting the presence of a distinct flow regime under the fast wave condition. Given the frequent occurrence of fast waves in the ocean environment, arising from the nonlinear interactions between slow waves or generated by local tsunami-induced waves propagating into the light-wind region (Åkervik & Vartdal 2019), the fast wave-induced momentum transfer requires further investigation and discussion.

The energy transfer between wind and surface waves has received considerable attention due to its critical role in determining wind–wave interaction and evolution. However, compared with the investigation of momentum transport above wind waves, there are fewer studies on kinetic energy transfer. Understanding the energy transfer involves analysing the conservation equation of kinetic energy, which includes mean, wave-induced and turbulent kinetic energy (TKE) components (Reynolds & Hussain 1972). The production terms in this conservation equation are essential for comprehending how energy is exchanged between different forms. Additionally, the terms associated with production, transport and dissipation play a significant role in balancing the kinetic energy. Despite efforts to study kinetic energy transfer above wind waves (Rutgersson & Sullivan 2005; Hara & Sullivan 2015; Yousefi, Veron & Buckley 2021), several questions still need to be answered.

Given the direct effect of TKE on energy transfer between wave-coherent and turbulent motions, previous studies paid much attention to the TKE variation, which was found to be relevant to the wave age (Shen et al. 2003; Sullivan et al. 2008). Since the turbulent flow over a slow wave is qualitatively similar to the stationary wavy wall turbulence, there is no doubt that the TKE is approximately consistent for both scenarios. Hara & Sullivan (2015) and Husain et al. (2019) noted that the TKE in a slow wave condition features a wave-coherent pattern. The intensified TKE observed on the leeward side of the wave is related to the strong shear layer behind the crest (Yang & Shen 2010). However, a higher wave age weakens the TKE peak and shifts it upstream (Shen et al. 2003; Buckley & Veron 2019), demonstrating the wave age dependency. Moreover, it is unsurprising that the TKE production is intensified behind the wave crest in a slow wave, the same as that of stationary wavy wall turbulence. Studies, through evaluating the TKE production profiles, suggest that the same trend can be found for both slow and fast wave regimes (Smedman et al. 1999): it is positive at all heights, reaches a peak near the wave surface and decreases vertically, agreeing basically with that of turbulent shear flow over a flat wall (Mansour, Kim & Moin 1988; Antonia et al. 1992). This indicates the energy transfer from wave-coherent shear to turbulence. It is known that the TKE production results from the interaction between different components of the work done by the turbulent stress on the velocity gradient. Therefore, further evaluation of these components should be conducted to physically reveal the TKE production mechanism.

The complexity of wind–wave energy transfer lies in the additional conservation equation of wave-induced kinetic energy (WKE), which describes the energy balance of wave-induced motion. In wind–wave interactions, the WKE production is vital in determining the energy transfer between wave-coherent and wave-induced motions. Hara & Sullivan (2015) pointed out that the WKE production cannot be ignored, although it is less than the TKE production, indicating the critical role of the WKE production in wind waves. It is also revealed that the WKE production is confined to the region near the surface wave (Rutgersson & Sullivan 2005; Hara & Sullivan 2015), mainly within the wave boundary layer (WBL). Yousefi et al. (2021) recently conducted experiments investigating the energy transfer above wind-generated waves and evaluated the WKE production in detail. They pointed out that the WKE production showed inverse features above and below the WBL, indicating different energy transfer mechanisms. It is inferred that the WKE or WKE production depends on wave age, as the wave-induced motion is strengthened with the increase of wave age. Due to the wave orbital velocity-induced motion, the WKE production is composed of the work done by wave-induced stress on the velocity gradient. However, how these components interact to determine the production of WKE under different wave regimes remains unclear.

Despite the production term, the WKE budget attaches a term depicting the wave–turbulence exchange, which also appears in the TKE budget. This term is characterized by the interaction between turbulent stress and wave-induced shear (Reynolds & Hussain 1972; Hara & Belcher 2004; Yousefi et al. 2021). Makin & Kudryavtsev (1999) pointed out that the TKE production consists of two sources: mean motion and wave-induced motion, with the latter presenting the wave–turbulence exchange. This indicates that the wave–turbulence exchange cannot be ignored when evaluating the macro-TKE production. Generally, the wave–turbulence exchange for a slow mechanically generated wave is positive above the height of $kz\approx 0.1$ (where $k$ is the wavenumber and $z$ is the vertical coordinate), leading to the energy transferred from wave to turbulence. However, according to the numerical results, Rutgersson & Sullivan (2005) found that, for a slow wave, the transport of wave–turbulence exchange is mainly in the opposite direction above the height of $kz\approx 0.1$, from turbulence to wave perturbation. They attributed the difference to the effect of the Reynolds number. The recent experimental results by Yousefi et al. (2021) suggest that the energy transfer between wave and turbulence shows an alternated positive–negative pattern along the wave; energy is transferred from turbulence to wave upwind of the wave crests and from wave to turbulence downwind. They also pointed out that a small amount of energy flows from turbulence to wave within the thin region above the wave. Notably, studies have put more effort into the wave–turbulence exchange for slow waves due to the complex dynamics. However, the fast wave-induced wave–turbulence exchanges still need further revelation.

It is worth noting that for wind–wave problems, a two-phase air–water flow model is the best way to replicate the wind–wave boundary layer flow due to the strong wind–wave coupling. However, these kinds of investigations mainly focus on wave breaking, dissipation and wave growth (Chen et al. 1999; Song & Sirviente 2004; Iafrati 2011; Deike, Popinet & Melville 2015; Wu & Deike 2021; Wu, Popinet & Deike 2022) while paying little attention to the turbulent flow dynamics. To our knowledge, the one-phase models were widely used in turbulent airflow over stationary wavy or prescribed moving wave boundaries (Sullivan et al. 2000; Kihara et al. 2007; Yang & Shen 2010; Druzhinin et al. 2012; Sullivan, McWilliams & Patton 2014; Sullivan et al. 2018; Cao & Shen 2021), providing deep insight into turbulent structures and momentum transfer above wind waves. Because of the strong coupling of the wind–wave problem, whether a one-way coupling method, such as the turbulent flow over moving boundaries, is suitable still needs to be clarified. In addition, the fast wave-induced momentum transfer mechanism and detailed energy transfer, including the TKE, WKE and wave–turbulence exchange between wind and waves under different wave regimes, need further revelation. Based on these motivations, the present paper simplifies wind waves as turbulent flow over a prescribing moving wave boundary concerning the momentum and energy transfers in the wind–wave interaction for different wave ages. The rest of the paper is structured as follows. Section 2 describes the mathematical formulation of the problem and numerical method. The momentum flux budget and Reynolds shear stress (RSS) production are presented in § 3. Section 4 shows the energy transfer in wind–wave interactions. The main conclusions are then summarized in § 5. For brevity and clarity, the derivations of the kinetic energy budgets and some verifications and discussions are included in the supplementary material available at https://doi.org/10.1017/flo.2024.22.

2. Mathematical formulation and numerical method

The present study simplifies the surface wave as a moving wave boundary. A small-amplitude surface wave propagating with phase velocity $c$ can be expressed as

(2.1)\begin{equation} \eta (x, z, t) = a\cos(k(x-ct)), \end{equation}

where $a$ is the wave amplitude and $k={2\pi }/{\lambda _x }$ denotes the wavenumber ($\lambda _x$ is the wavelength in the direction of propagation). The turbulent flow above the propagating wave can be described through the three-dimensional incompressible Navier–Stokes equations

(2.2)$$\begin{gather} \frac{\partial {{u}_{i}}}{\partial {{x}_{i}}} = 0, \end{gather}$$

(2.3)$$\begin{gather}\frac{\partial {{u}_{i}}}{\partial t} + {{u}_{j}}\frac{\partial {{u}_{i}}}{\partial {{x}_{j}}} = {-}\frac{1}{\rho }\frac{\partial p}{\partial {{x}_{i}}} + \nu \frac{{{\partial }^{2}}{{u}_{i}}}{\partial {{x}_{j}}\partial {{x}_{j}}} + \varPi {{\delta }_{1i}}, \end{gather}$$

where $(x_1,x_2,x_3)=(x,y,z)$ (with tensor notation $i$ the $i$th component), respectively, denote the streamwise, spanwise and vertical coordinates, $(u_1,u_2,u_3)=(u,v,w)$ represent the corresponding velocity components, and repeated indices are implicitly summed over. Here $p$ is the pressure, $\varPi$ is the external force transformed as a pressure gradient driving the flow, ${{\delta }_{1i}}$ is the Kronecker delta with ${{\delta }_{11}}=1$ and ${{\delta }_{12}}={{\delta }_{13}}=0$, $\nu$ is the kinematic viscosity and $\rho$ is the fluid density.

2.1 Mesh motion equation

The mesh motion equations are coupled to the conservation equations of fluid motion to simulate the moving wavy boundary. It is assumed that the disturbance of the propagating wave on the whole field decreases with the increase of the distance from the moving boundary, and the disturbance at infinity is zero. Therefore, for the entire computational domain, the mesh motion can be described by the diffusion equation governing the grid point displacement,

(2.4)\begin{equation} K\frac{\partial^2D_{w,j}}{\partial x_i\partial x_i} = 0,\quad\text{for } j = 1,2,3, \end{equation}

where $({{D}_{w,1}},\,{{D}_{w,2}},\,{{D}_{w,3}})$ represent the displacement components for the $w$th grid point in the streamwise, spanwise and vertical directions, respectively. Here, $K$ denotes the diffusive coefficient related to the position of the grid point that can be defined by $K={{l}^{-w}}$, where $l$ is the average distance from the wall to every grid point, namely the inverse distance.

2.2 The Navier–Stokes equations in a moving frame

The grid point displacement at the boundary can be determined according to (2.1). Hence, the displacement of each grid point can be calculated by (2.4). The displacement vector for the $w$th grid point is $D_{w,i}^{t+\Delta t}-D_{w,i}^{t}$ from $t$ to $t+\Delta t$, and the velocity vector of the moving grid point is thus ${{u}_{w,i}}=({D_{w,i}^{t+\Delta t}-D_{w,i}^{t}})/{\Delta t}$. The fluid motion coupled with the mesh motion can be described by the Navier–Stokes equations in a moving frame:

(2.5)$$\begin{gather} \frac{\partial ({{u}_{i}}-{{u}_{w,i}})}{\partial {{x}_{i}}} = 0, \end{gather}$$

(2.6)$$\begin{gather}\left.\frac{\partial {{u}_{i}}}{\partial t}\right|_{\varsigma } + ( {{u}_{j}}-{{u}_{w,j}} )\frac{\partial {{u}_{i}}}{\partial {{x}_{j}}} = {-}\frac{1}{\rho }\frac{\partial p}{\partial {{x}_{i}}} + \nu \frac{{{\partial }^{2}}{{u}_{i}}}{\partial {{x}_{j}}\partial {{x}_{j}}} + \varPi {{\delta }_{1i}}. \end{gather}$$

Compared with the equations in the stationary coordinate frame, the mesh motion would introduce the time change rate of the velocity in the new coordinate frame $\varsigma$, which is independent of the Lagrangian and Eulerian coordinate systems. The mesh motion also introduces the velocity vector of the moving coordinate frame ${\boldsymbol {u}_{w}}=( {{u}_{w,1}},\,{{u}_{w,2}},\,{{u}_{w,3}} )$ in the mass conservation equation and the nonlinear terms of the momentum conservation equations.

The integral conservation equation of a general tensorial property $\psi$ on a control volume $V$ bounded by a closed surface $S$ with arbitrary motion is

(2.7)\begin{equation} \frac{{\rm d}}{{\rm d}t}\int_{V} \rho \psi\,{\rm d}V + \oint_{S}{\rho \boldsymbol{n}\boldsymbol{\cdot} ( \boldsymbol{u}-{{\boldsymbol{u}}_{w}} )}\psi \,{\rm d}S = {-}\oint_{S}{\rho \boldsymbol{n}\boldsymbol{\cdot} {{\boldsymbol{q}}_{\psi }}\,{\rm d}S} + \int_{V}{{{S}_{\psi }}\,{\rm d}V}, \end{equation}

where $\boldsymbol {u}$ is the fluid velocity, ${{\boldsymbol {u}}_{w}}$ is the velocity vector of the control volume, ${{\boldsymbol {q}}_{\psi }}$ denotes the volume source, ${{S}_{\psi }}$ is the surface source and $\boldsymbol {n}$ is the unit normal vector pointing outward on the boundary surface of a grid (Jasak & Tuković 2006). Replacing $\psi$ with the fluid velocity vector $\boldsymbol {u}$ can lead to the integral conservation momentum equations under the moving coordinate frame,

(2.8)\begin{equation} \frac{{\rm d}}{{\rm d}t}\int_{V}{\rho \boldsymbol{u}\,{\rm d}V} + \oint_{S}{\rho \boldsymbol{n}\boldsymbol{\cdot} ( \boldsymbol{u}-{{\boldsymbol{u}}_{w}} )}\boldsymbol{u}\,{\rm d}S = {-}\oint_{S}{\rho \boldsymbol{n}\boldsymbol{\cdot} {{\boldsymbol{q}}_{\boldsymbol{u}}}\,{\rm d}S} + \int_{V}{{{S}_{\boldsymbol{u}}}\,{\rm d}V}. \end{equation}

The second-order finite volume discretization of (2.8) for a structural grid, achieved using the midpoint rule, can convert the surface integral into the sum of area parts. After discretization, the momentum equations can be expressed as

(2.9)\begin{equation} \frac{{{( {{\rho }_{P}}{{\boldsymbol{u}}_{P}}{{V}_{P}} )}^{t + \Delta t}}-{{( {{\rho }_{P}}{{\boldsymbol{u}}_{P}}{{V}_{P}} )}^{t}}}{\Delta t} + \sum_{f}{{{\rho }_{f}}( F-{{F}_{S}} ){{\boldsymbol{u}}_{f}}} = {-}\sum_{f}{{{\boldsymbol{S}}_{f}}\boldsymbol{\cdot} \rho {{\boldsymbol{q}}_{\boldsymbol{u}}}} + {{S}_{\boldsymbol{u}}}{{V}_{P}}, \end{equation}

where the subscript $P$ denotes the grid point, $f$ the face value, ${{V}_{P}}$ the grid volume, $F={{\boldsymbol {S}}_{f}}\boldsymbol {\cdot } {{\boldsymbol {u}}_{f}}$ the fluid flux with ${{\boldsymbol {S}}_{f}}=\boldsymbol {n}{{S}_{f}}$, and the mesh motion flux is ${{F}_{S}}$. As described in (2.9), the introduced time change rate for the time term is converted to the mesh motion flux, which should meet the need for space conservation,

(2.10)\begin{equation} \frac{{\rm d}}{{\rm d}t}\int_{V}{\,{\rm d}V}-\oint_{S}{\boldsymbol{n}\boldsymbol{\cdot} {{\boldsymbol{u}}_{w}}}\,{\rm d}S = 0. \end{equation}

Equation (2.10) should be satisfied for every time step, with its discretizing form expressed as

(2.11)\begin{equation} \frac{V_{P}^{t + \Delta t}-V_{P}^{t}}{\Delta t}-\sum_{f}{{{F}_{S}}} = 0. \end{equation}

It is noted that the mesh motion flux ${{F}_{S}}$ is determined by the volume swept by the moving surface $f$ during the current time step, not by the mesh velocity ${{\boldsymbol {u}}_{w}}$.

The present paper uses the finite volume method to solve the governing equations. The OpenFOAM (open field operation and manipulation) solver is used to accomplish the simulation. It is noted that we here directly solve the full Navier–Stokes equations without a turbulent model, which can be regarded as an approximation of direct numerical simulation (DNS), termed the quasi-DNS. The DNS requires fine mesh, resolving the Kolmogorov scale. However, quasi-DNS loosens this constraint, as verified by Komen et al. (2017) (figures 3 and 5 in their paper). Therefore, to accomplish the simulation accurately, we use a relatively fine mesh (discussed in the supplementary material for grid convergence). In the current solver, the convection terms in the momentum equations are spatially discretized by the second-order upwind difference scheme, and the second-order backward implicit scheme is adopted for time stepping. The finite volume method based on the standard pressure velocity coupling method and PIMPLE algorithm is a variant of the PISO (pressure-implicit with splitting of operators) method to solve the governing equations. The PIMPLE algorithm combines the SIMPLE (semi-implicit method for pressure-linked equations) strategy and the PISO algorithm and regards each time step as a steady-state flow (Gatin et al. 2017; Gimenez et al. 2019). The standard PISO algorithm is used for the last step when the solution is obtained according to the steady-state algorithm to a certain extent.

2.3 Simulation configuration

The turbulent flow over a propagating wave is simplified as channel flow with a translating wavy bottom, driven by the external force transformed as a pressure gradient to ensure the fixed bulk velocity ${{U}_{0}}={\iint {U\,{\rm d}x\,{\rm d}z}}/{\iint {{\rm d}x\,{\rm d}z}}$, where $U$ is the instantaneous velocity, which is the same as the numerical study by Åkervik & Vartdal (2019). The Reynolds number based on the bulk velocity and wavelength is $Re={{{U}_{0}}{{\lambda }_{x}}}/{\nu }\approx 7300$ to provide fully developed turbulence. The corresponding wave friction Reynolds numbers ${Re}_{\tau }={{{u}_{\tau }}{{\lambda }_{x}}}/{\nu }$ are shown in table 1. In the wind–wave turbulent boundary layer flow, the wave slope and wave age determine the dynamics of wind–wave interactions. The former is defined by the ratio of amplitude $a$ to the wavelength ${{\lambda }_{x}}$ (or the product of amplitude $a$ and wavenumber $k$). It is noted that the present paper considers the Airy waves with the phase speed $c$, as described in (2.1). The Airy wave is the linear water wave solution and requires a small wave slope (Yang & Shen 2010). The wave slope in the present study is thus $ak=0.1$, 0.13, 0.15, 0.20, 0.25. Three regimes of wind–wave conditions are considered: slow (${c}/{{{u}_{\tau }}}=2.01$, 3.28, 3.43, 3.69); intermediate (${c}/{{{u}_{\tau }}}=16.36$, 16.77, 17.25); and fast (${c}/{{{u}_{\tau }}}=35.80$, 42.83, 58.93, 63.65) waves. The size of the computational domain is $( x,y,z )=( 2{{\lambda }_{x}},{{\lambda }_{x}},0.5{{\lambda }_{x}})$, which is verified to be large enough (see the supplementary material). The grid points are evenly spaced in both streamwise and spanwise directions. In the vertical direction, the grid points are clustered at the wavy wall boundary through exponential transformation to enhance the accuracy for the boundary layer, in which the first layer of the grid meets the need of $\Delta \zeta _{w}^{+}<1$, where + denotes the normalization scaled by viscous unit, namely $\Delta \zeta _{w}^{+}={\Delta {{\zeta }_{w}}{{u}_{\tau }}}/{\nu }$. Here, we used wave-following curvilinear coordinate $(\xi,\zeta )=(x,{(z-\eta )}/{H})$ to show the dimensionless grid scales, where $\eta$ is the wave elevation and $H$ denotes the height of the physical domain (Yang & Shen 2017; Cao, Deng & Shen 2020). The high grid resolution (HGR) and super high grid resolution (SHGR) cases are used for grid convergence, and cases S4 and S5 are used to compare with the numerical results by Yang & Shen (2010, 2017), discussed in the supplementary material. The total number of grid points for case SHGR is ${{N}_{x}}\times {{N}_{y}}\times {{N}_{z}}=401\times 201\times 151$, while the other cases use ${{N}_{x}}\times {{N}_{y}}\times {{N}_{z}}=251\times 126\times 76$, with the dimensionless grid scale being shown in table 1. It is seen that the grid scale reaches a quasi-DNS (Komen et al. 2017) but cannot reach the Kolmogorov scale. Still, the grid convergence in the supplementary material demonstrates the suitability of the chosen grid scale. The periodic conditions are applied along the streamwise and spanwise directions. The propagating wave boundary is applied to the bottom wall, while the upper wall is set to be a slip wall. After the full development of the turbulent flow, the statistical averaging begins.

Table 1. Parameter settings. The HGR and SHGR cases are used for the verification of grid convergence. The present study considers five wave slopes: $ak$ = 0.1, 0.13, 0.15, 0.2, 0.25. The cases are divided into slow wave (S1–S5), intermediate wave (I1–I3) and fast wave (F1–F4).

2.4 Phase-averaging decomposition

To separate the effects of wave motion and turbulent motion, we use the phase-averaging method to obtain the disturbance feature caused by the waves, which had been widely used in previous investigations of wind waves (Sullivan et al. 2000; Kihara et al. 2007; Yang & Shen 2010, 2017; Yousefi et al. 2020a, 2021; Wu et al. 2022). The instantaneous quantity can be decomposed into phase-averaging and fluctuating components, and the phase-averaging quantity (also called the wave-coherent quantity) can further be decomposed into ensemble-averaged and wave-induced quantities. The wave-induced quantity represents the disturbance of waves to the above macroscopic flow. Hence, the instantaneous quantity can be decomposed into three parts (triple decomposition):

(2.12)$$\begin{gather} {{f}_{i}}( x,z,t ) = \langle\,{{f}_{i}}( \varphi ,z )\rangle + {{{f}'_i}}( x,z,t ) = \overline{{{f}_{i}}}( z ) + \widetilde{{{f}_{i}}}( \varphi ,z ) + {{f}'_{i}}( x,z,t ), \end{gather}$$

(2.13)$$\begin{gather}\langle\, {{f}_{i}}( \varphi ,z )\rangle = \underset{N\to \infty }{\mathop{\lim }}\,\frac{1}{N}\sum_{n = 0}^{N-1}{{{f}_{i}}( \varphi + 2n\pi ,z )}, \end{gather}$$

(2.14)$$\begin{gather}\overline{{{f}_{i}}}( z ) = \underset{M\to \infty }{\mathop{\lim }}\,\frac{1}{M}\sum_{m = 0}^{M-1}\left\langle\, {{f}_{i}}\left( \varphi + \frac{2\pi m}{M},z \right)\right\rangle. \end{gather}$$

Here $\widetilde {{{f}_{i}}}( \varphi,z )$ is the wave-induced quantity, determined by the difference between the phase-averaging quantity $\langle\, {{f}_{i}}( \varphi,z )\rangle$ and the ensemble-averaged quantity $\overline {{{f}_{i}}}( z )$, and ${{f}'_{i}}( x,z,t )$ represents the turbulent fluctuation. It is worth noting that the present study only considers the single wave. The triple decomposition is a suitable way to separate the flow fields. In the real marine environment, the broad-banded wave fields are common, which means the failure of this decomposition. Hristov, Friehe & Miller (1998) pointed out that the waves in the open ocean are continuously spread throughout spectral scales. Taking a period of the most energetic (peak) mode in the wave spectrum as the characterized period and conducting phase-averaging leads to a strong attenuation due to the destructive interference of multiple modes. They proposed a novel approach, employing an eikonal-like representation of the wave field based on the concept of an analytic signal and the Hilbert transform to identify the wave-coherent component. For monochromatic signals of wave profiles, the results filtered by this method are equivalent to those obtained by phase-averaging.

2.5 Verification of the numerical model

To verify the numerical model, figure 1 shows the comparisons of momentum statistics with the experimental results by Yousefi et al. (2020a, b) and numerical results by Yang & Shen (2010, 2017). All quantities are scaled by the friction velocity. We have used the friction velocity determined by extrapolating the ensemble-averaged RSS from the outer layer based on the physical essence of the friction velocity. The downward shifted logarithmic and linear laws compared with the flat-wall boundary layer flow (FWBL) in figure 1(a) agree well with the experimental and other simulation results. In the supplementary material, we also compare the results of channel flow with the direct numerical simulation results by Lee & Moser (2015, 2018). Figure 1(b) shows the balance of momentum flux for case S1. The RSS and wave-induced shear stress (WSS) trends agree qualitatively with Åkervik & Vartdal (2019) (pressure-driven turbulent flow), although the wave friction Reynolds number ($Re_{\tau }=557$) for case S1 is different from that in Åkervik & Vartdal (2019) ($Re_{\tau }=200$, 260, 395, 950). It is noted that the Reynolds number is an essential parameter for the scenarios of wind–wave interactions. Åkervik & Vartdal (2019) provides insights into the friction Reynolds number effect, which covers a range of wave ages and Reynolds numbers, that the budget of momentum flux shows similar pattern for low- and high-friction Reynolds numbers for different wave ages. In fact, the wave age can be written as a ratio of Reynolds number, namely ${c}/{{{u}_{\tau }}}=2\pi ({c}/{k\nu })({\nu }/{{{u}_{\tau }}\lambda _{x}})=2\pi ({Re_{w}}/{Re_{\tau }})$, where $Re_{w}$ is the wave Reynolds number. Since the present results are normalized by the friction velocity, we do not pay much attention to the friction Reynolds number effect. It is also found that the viscous shear stress (VSS) agrees well with Yousefi et al. (2020b). Figure 1(c) indicates a slight discrepancy compared with experimental results, which is attributed to the different coordinate frame, as reported in figure 13 by Yousefi et al. (2020a). Figure 1(d) shows consistency with the numerical results by Yang & Shen (2010). It is noted that Yang & Shen (2010) considered the stress-driven turbulent Couette flow over a progressive surface wave train, which differs from the pressure-driven turbulent flow in the current study. Therefore, the RSS is close to the wall stress of the upper plane, whereas the present RSS is nearly reduced to zero. Still, the RSS variation in the near wave region is approximately consistent for both cases. Furthermore, we plot the RSS and VSS along the wave propagation direction at the same height as Yang & Shen (2010) and Yousefi et al. (2020a) (the RSS locates at the height ${\zeta }/{{{\lambda }_{x}}\approx 0.1}$, and the VSS is obtained by averaging stress measurements between 284 and 664 $\mathrm {\mu }$m above the air–water interface corresponding to the region of ${{\zeta }^{+}}\approx 3.38$–7.5 in the present study), as shown in figure 1(e,f). We see both RSS and VSS basically agree with the previous investigations. It is noted that we use a common method to determine the friction velocity adopted in previous studies (Hamed et al. 2015; Bopp 2018), that the total wind stress can be determined by extrapolating turbulent stress from the outer layer since it matches the sum of turbulent, wave-coherent and viscous stress in that region. To address this, we simulated additional slow wave cases to highlight the wind stress partition. Figure 2 shows the form drag and viscous stress normalized by $\rho u_{\tau }^2$, where $u_{\tau }$ is determined by this method. For comparison, the data from many previous studies are also plotted (Banner 1990; Banner & Peirson 1998; Sullivan et al. 2000; Caulliez et al. 2008; Peirson & Garcia 2008; Savelyev 2009; Grare et al. 2013b; Peirson et al. 2014; Bopp 2018; Sullivan et al. 2018; Buckley et al. 2020; Funke et al. 2021). Figure 2 suggests that the mean form drags fall within the region overlapping with findings from Peirson et al. (2014) and Funke et al. (2021). The viscous stress also falls into the region very close to that observed in previous studies. Moreover, it approaches the dashed line determined by subtracting form drag from the total wind stress. Therefore, the extrapolating RSS from the outer layer to the surface waves can be a viable approach for calculating total wind stress. By evaluating these turbulent high-order statistics, form drag and viscous stress, it is inferred that the present numerical model can provide rational and reliable data.

Figure 1. Comparisons of the momentum statistics between the present study and other numerical and experimental results: EXP, experimental results by Yousefi et al. (2020a) and Yousefi, Veron & Buckley (2020b) for $ak=0.13$ and ${c}/{{{u}_{\tau }}}=3.69$; DNS1, numerical results by Yang & Shen (2017) for $ak=0.25$ and ${c}/{{{u}_{\tau }}}=2$; DNS2, numerical results by Yang & Shen (2010) for $ak=0.1$ and ${c}/{{{u}_{\tau }}}=2$. (a) Velocity profiles (the law of flat wall boundary layer flow is also shown by the black doubled dotted dashed line). (b) Vertical profiles of the shear stress (momentum flux) components (the dashed line denotes the extrapolation of RSS from the outer layer). Comparisons of RSS and WSS with (c) EXP and (d) DNS results. The discrepancy with EXP is ascribed to the different coordinate frames (it can be seen in figure 13 of Yousefi et al. (2020a)). (e) The comparison of RSS with DNS2 along the wave propagation direction, located at the height ${\zeta }/{{{\lambda }_{x}}\approx 0.1}$. (f) The comparison of VSS with EXP along the wave propagation direction.

Figure 2. (a) Mean form drag and (b) viscous stress. The data from previous studies are also plotted (Banner 1990; Banner & Peirson 1998; Sullivan et al. 2000; Caulliez, Makin & Kudryavtsev 2008; Peirson & Garcia 2008; Savelyev 2009; Grare et al. 2013b; Peirson, Walker & Banner 2014; Bopp 2018; Sullivan et al. 2018; Buckley, Veron & Yousefi 2020; Funke et al. 2021).

3. Momentum transfer and RSS production

Figure 3 shows the budgets of momentum flux for cases S1, I1 and F1, along with the extrapolation of RSS profiles from the outer layer. As the wave age increases, the RSS decreases, whereas the WSS transitions from negative to slightly positive values and eventually becomes strongly positive. This transition is highly related to the wave vertical motion (see the supplementary material). The turbulent flow over slow waves shares in common with the stationary two-dimensional wavy wall turbulence. Using the dispersive shear stress filtered by double-averaging (spatial and temporal averaging) in rough wall turbulence for reference, we can see some similarities, that the negative WSS in figure 3(a) suggests a form effect, approximately equivalent to the effect of dispersive shear stress in studies of rough wall turbulence. It is noted that Dey et al. (2020) considered the flow over two-dimensional dunes, which is relatively like the present wavy surface. We observe a same trend of dispersive shear stress as the present slow WSS. However, as the wave propagates quickly, the wave motion can provide enough perturbation to disturb the upper turbulent flow. The WSS becomes positive to achieve the momentum transfer from wave to wind. Therefore, we infer that from slow to fast waves, the transition from form effect to motion effect dominates the variation of WSS. The VSS also intensifies in cases with high wave age and dominates within the region of ${{\zeta }^{+}}<10$. The variation of momentum flux demonstrates that fast waves degenerate RSS, with strengthened WSS and VSS balancing this deficit.

Figure 3. Budgets of the shear stress, including the RSS, WSS and VSS (the dashed lines denote the extrapolation of the RSS profiles from the outer layer). Here (a) S1 with $ak=0.13$ and $c/u_{\tau }=3.69$; (b) I1 with $ak=0.13$ and $c/u_{\tau }=17.25$; (c) F1 with $ak=0.13$ and $c/u_{\tau }=35.80$.

Furthermore, we discuss how the RSS is produced. The RSS production is written as (see the supplementary material)

(3.1)\begin{equation} {{\mathcal{P}}_{13}} = \underbrace{-\langle {u}'{u}'\rangle \frac{\partial \langle w\rangle }{\partial x}}_{{{R}_{1131}}}\underbrace{-\langle {u}'{w}'\rangle \frac{\partial \langle w\rangle }{\partial z}}_{{{R}_{1333}}}\underbrace{-\langle {u}'{w}'\rangle \frac{\partial \langle u\rangle }{\partial x}}_{{{R}_{1311}}}\underbrace{-\langle {w}'{w}'\rangle \frac{\partial \langle u\rangle }{\partial z}}_{{{R}_{3313}}}, \end{equation}

where $R_{1131}$ denotes the streamwise Reynolds normal stress (SRNS) performing work on the streamwise gradient of vertical velocity, $R_{1333} (R_{1311})$ denotes the work done by RSS on the vertical gradient of vertical velocity (streamwise gradient of streamwise velocity), and $R_{3313}$ denotes the vertical Reynolds normal stress (VRNS) performing work on the vertical gradient of streamwise velocity.

Figure 4 shows the contour of $\mathcal {P}_{13}^{+}$ and its components scaled by ${u_{\tau }^{4}}/{\nu }$. It shows that as a primary factor, $R_{3313}^{+}$ contributes to the RSS production in case S1 in figure 4(a). In other words, the VRNS interacts with the vertical gradient of the wave-coherent streamwise velocity to determine the RSS production in a slow wave condition. For case I1, as shown in figure 4(b), both $R_{1333}^{+}$ and $R_{3313}^{+}$ negatively contribute to the RSS production near the crest, whereas $R_{1131}^{+}$ intensifies above the trough. As shown in figure 4(c), fast waves further amplify this pattern.

Figure 4. The RSS production and its components for cases (a) S1 with $ak=0.13$ and $c/u_{\tau }=3.69$; (b) I1 with $ak=0.13$ and $c/u_{\tau }=17.25$; (c) F1 with $ak=0.13$ and $c/u_{\tau }=35.80$. The red dashed line is the critical height. Here, $R_{1131}^{+}$ is the contribution by SRNS performing work on the streamwise gradient of vertical velocity, $R_{1333}^{+}$ ($R_{1311}^{+}$) is the contribution by RSS performing work on the vertical gradient of vertical velocity (streamwise gradient of streamwise velocity) and $R_{3313}^{+}$ is the contribution by VRNS performing work on the vertical gradient of streamwise velocity.

Figure 5 shows the ensemble-averaged profiles of $\bar {\mathcal {P}}_{13}^{+}$ and its components non-dimensionalized by ${u_{\tau }^{4}}/{\nu }$. Figure 5(a) indicates that the terms related to the SRNS and RSS can be ignored for slow waves. This is approximately equivalent to two-dimensional shear turbulence with a positive mean strain rate since the RSS is produced via the interaction between VRNS and the vertical gradient of the streamwise velocity. With the increase of wave age, as shown in figure 5(b), $\bar {\mathcal {P}}_{13}^{+}$ is reduced and approximately governed by the interaction between $\bar {R}_{1131}^{+}$, $\bar {R}_{1333}^{+}$ and $\bar {R}_{3313}^{+}$, with the strength of $\bar {R}_{1131}^{+}$ (favourable to the RSS production) lower than the sum of $\bar {R}_{1333}^{+}$ and $\bar {R}_{3313}^{+}$ (both of which are against the RSS production). However, a fast wave, as shown in figure 5(c), enhances the component $\bar {R}_{1131}^{+}$, thus leading to the rapid transformation of $\bar {\mathcal {P}}_{13}^{+}$ along the vertical direction. It is noted that the contour results in figure 4 do not clearly reflect the contribution of these components especially for the intermediate and fast waves. We attribute this to the ensemble-averaging conduction, which would eliminate some contributions of positive $R_{1131}^+$. Therefore, the contribution of $\bar {R}_{1131}^{+}$ is not much higher than that of $\bar {R}_{1333}^{+}$ and $\bar {R}_{3313}^{+}$.

Figure 5. The profiles of RSS production and its components for cases (a) S1 with $ak=0.13$ and $c/u_{\tau }=3.69$; (b) I1 with $ak=0.13$ and $c/u_{\tau }=17.25$; (c) F1 with $ak=0.13$ and $c/u_{\tau }=35.80$. The solid line denotes the total RSS production, the dashed line is the contribution of $\bar {R}_{1131}^{+}$, the dotted line is the contribution of $\bar {R}_{1333}^{+}$, the dotted dashed line is the contribution of $\bar {R}_{1311}^{+}$, and the double-dotted dashed line is the contribution of $\bar {R}_{3313}^{+}$.

It can be inferred that the RSS production for a slow wave exhibits similarities to two-dimensional shear turbulence with a positive mean strain rate. However, as the wave age increases, the RSS is produced in a complex manner between SRNS, VRNS and RSS-related production terms. The intermediate wave age can be regarded as a transition because $\bar {R}_{1131}^{+}$ cannot counteract $\bar {R}_{1333}^{+}$.

4. Energy transfer in wind–wave interactions

4.1 Energy transfer between wave-coherent, and turbulent and wave-induced motions (TKE or WKE production)

This section discusses the energy transfer between phase-averaged (wave-coherent), wave-induced and turbulent motions. We here use the two-dimensional statistical field; consequently, $T_{t}$ (denotes the energy transfer between wave-coherent and turbulent motions) and $T_{w}$ (denotes the energy transfer between wave-coherent and wave-induced motions) can be expanded as (see the supplementary material)

(4.1)$$\begin{gather} {{T}_{t}} = \underbrace{-\langle {u}'{u}'\rangle \frac{\partial \langle u\rangle }{\partial x}}_{{{T}_{t,11}}}\underbrace{-\langle {u}'{w}'\rangle \frac{\partial \langle u\rangle }{\partial z}-\langle {u}'{w}'\rangle \frac{\partial \langle w\rangle }{\partial x}}_{{{T}_{t,13}}}\underbrace{-\langle {w}'{w}'\rangle \frac{\partial \langle w\rangle }{\partial z}}_{{{T}_{t,33}}}, \end{gather}$$

(4.2)$$\begin{gather}{{T}_{w}} = \underbrace{-\tilde{u}\tilde{u}\frac{\partial \langle u\rangle }{\partial x}}_{{{T}_{w,11}}}\underbrace{-\tilde{u}\tilde{w}\frac{\partial \langle u\rangle }{\partial z}-\tilde{u}\tilde{w}\frac{\partial \langle w\rangle }{\partial x}}_{{{T}_{w,13}}}\underbrace{-\tilde{w}\tilde{w}\frac{\partial \langle w\rangle }{\partial x}}_{{{T}_{w,33}}}, \end{gather}$$

where ${{T}_{t,11}}$ (${{T}_{t,13}}$ or ${{T}_{t,33}}$) is TKE production through SRNS (RSS or VRNS) performing work on the velocity gradient, and ${{T}_{w,11}}$ (${{T}_{w,13}}$ or ${{T}_{w,33}}$) is WKE production through wave-induced streamwise normal stress (WSS or wave-induced vertical normal stress) performing work on velocity gradient.

Figure 6 shows the total and components of TKE production caused by the interaction between turbulent stresses and phase-averaged shear (also known as the energy transfer between wave-coherent and turbulent motions), with all components normalized by ${u_{\tau }^{4}}/{\nu }$. A slow wave in figure 6(a) generates locally enhanced turbulent production on the leeward side. Conversely, slightly negative production occurs above the windward side, where acceleration induces a favourable pressure gradient, which agrees well with Yousefi et al. (2021). However, an intermediate wave in figure 6(e) sharply decreases this production on the leeward side while intensifying it on the windward side. A fast wave in figure 6(i) further strengthens the turbulent production on the windward side but decreases it into a negative value on the leeward side, which denotes turbulent consumption, where energy is transferred from turbulence to wave-coherent motion.

Figure 6. Energy transfer between wave-coherent and turbulent motions. The total and components of TKE production ($T_{t}^{+},T_{t,11}^{+},T_{t,13}^{+},T_{t,33}^{+}$) are plotted for (a–d) slow wave (S1 with $ak=0.13$ and $c/u_{\tau }=3.69$), (e–h) intermediate wave (I1 with $ak=0.13$ and $c/u_{\tau }=17.25$) and (i–l) fast wave (F1 with $ak=0.13$ and $c/u_{\tau }=35.80$). The red dashed line is the critical height.

For the wave ages considered in the present study, $T_{t,11}^{+}$ and $T_{t,13}^{+}$ dominate the total production. Specifically, a negative (positive) $T_{t,11}^{+}$ on the leeward (windward) side and a positive (negative) $T_{t,13}^{+}$ on the leeward (windward) side can be observed across different wave regimes. It is noted that Yousefi et al. (2021), by orthogonal curvilinear coordinate frame, showed that $T_{t,13}^{+}$ is always positive. But our simulation results agree with those of Hudson et al. (1996), Yang & Shen (2010, 2017) and Buckley & Veron (2019) under the Cartesian coordinate frame. By adjusting the relative magnitude of $T_{t,11}^{+}$ and $T_{t,13}^{+}$, the waves determine whether turbulence is generated or consumed. In figure 6(a–c), $T_{t,11}^{+}$ competes with $T_{t,13}^{+}$, with $T_{t,11}^{+}< T_{t,13}^{+}$ ($T_{t,11}^{+}>T_{t,13}^{+}$) prevailing on the leeward (windward) side in a slow wave. But when the wave becomes fast, as shown in figure 6(j,k), $T_{t,11}^{+}$ is always stronger than $T_{t,13}^{+}$ on both sides. In other words, $T_{t,13}^{+}$ negatively contributes to $T_{t}^{+}$, with normal stress playing a vital role in the energy transfer between wave-coherent and turbulent motions. The intermediate wave in figure 6(f,g) can be regarded as the transitional condition, where the competition starts to change. Overall, a fast wave enhances the streamwise Reynolds normal stress, the windward side's negative RSS and the velocity gradient, as shown in the supplementary material, to produce additional energy transfer between wave-coherent and turbulent motions.

To further investigate the energy transfer between wave-coherent and wave-induced motions, we present the contours of total and components of WKE production (non-dimensionalized by ${u_{\tau }^{4}}/{\nu }$) in figure 7. The positive (negative) total WKE production $T_{w}^{+}$ denotes that the wave-induced motion produces (disrupts) kinetic energy. Figure 7(a) shows that wave-induced motion generates (disrupts) kinetic energy on the windward (leeward) side for a slow wave. However, the rising wave age in figure 7(e,i) results in a relatively symmetrical distribution and strengthened WKE production.

Figure 7. Energy transfer between wave-coherent and wave-induced motions. The total and components of WKE production ($T_{w}^{+},T_{w,11}^{+},T_{w,13}^{+},T_{w,33}^{+}$) are plotted for (a–d) slow wave (S1 with $ak=0.13$ and $c/u_{\tau }=3.69$), (e–h) intermediate wave (I1 with $ak=0.13$ and $c/u_{\tau }=17.25$), and (i–l) fast wave (F1 with $ak=0.13$ and $c/u_{\tau }=35.80$). The red dashed line is the critical height.

The components of WKE production in figure 7 reveal an intriguing energy production mechanism. A slow wave leads to a typical interaction between $T_{w,11}^{+}$ and $T_{w,13}^{+}$, as shown in figure 7(b,c), where $T_{w}^{+}$ is determined through the near wave's WSS and the upper wave-induced streamwise normal stress, with the rare $T_{w,33}^{+}$ contribution in figure 7(d). This agrees with Yousefi et al. (2021), who found that $T_{w,11}^{+}$ and $T_{w,13}^{+}$ contribute to the total wave production for high winds (corresponding to slow waves). Nevertheless, a fast wave transforms another interaction mechanism shown in figure 7(j–l): both $T_{w,13}^{+}$ and $T_{w,33}^{+}$ contribute to WKE production. Similarly, the intermediate wave in figure 7(f–h) is treated as the transitional condition, where this interaction mechanism starts to change.

4.2 Energy transfer between turbulence and waves (wave–turbulence exchange)

As described in the supplementary material, ${{W}_{t}}$ represents the energy transfer between wave and turbulence (Reynolds & Hussain 1972; Hara & Belcher 2004; Yousefi et al. 2021), which means waves produce or disrupt the TKE, called wave–turbulence exchange. While ${{W}_{w}}$ denotes the WKE transportation, indicating the wave motion generated by the wave-induced shear. As emphasized by Yousefi et al. (2021), ${{\bar {W}}_{w}}$ redistributes wave-induced motions. It can be incorporated into the transport term of the kinetic energy density of the wave-induced flow. Some previous studies suggest that ${{\bar {W}}_{w}}$ is negligible compared with other transport terms (Makin & Kudryavtsev 1999; Hara & Belcher 2004; Hara & Sullivan 2015). Therefore, we only show ${{W}_{t}}$ here and the WKE transportation ${{W}_{w}}$ is discussed in the supplementary material. Following § 4.1 and expanding $W_t$ as

(4.3)\begin{equation} {{W}_{t}} = \underbrace{-\langle {u}'{u}'\rangle \frac{\partial \tilde{u}}{\partial x}}_{{{W}_{t,11}}}\underbrace{-\langle {u}'{w}'\rangle \frac{\partial \tilde{u}}{\partial z}-\langle {u}'{w}'\rangle \frac{\partial \tilde{w}}{\partial x}}_{{{W}_{t,13}}}\underbrace{-\langle {w}'{w}'\rangle \frac{\partial \tilde{w}}{\partial z}}_{{{W}_{t,33}}}, \end{equation}

where ${{W}_{t,11}}$ (${{W}_{t,13}}$ or ${{W}_{t,33}}$) is wave–turbulence exchange through SRNS (RSS or VRNS) performing work on the wave-induced velocity gradient.

Figure 8 shows $W_{t}^{+}$ and its components rescaled by ${u_{\tau }^{4}}/{\nu }$. A positive value represents that the energy is transferred from waves to turbulence, whereas a negative value means the other way around. It is typical for all wave regimes that waves inject (extract) energy into (from) turbulence on the windward (leeward) side, as shown in figure 8(a,e,i). Notably, the enhanced positive $W_{t}^{+}$ on the windward side for a slow wave extends downstream, resulting from the combination of turbulent stress and wave-induced strain rate. Figure 8(a–d) present that a slow wave determines $W_{t}^{+}$ primarily through $W_{t,11}^{+}$ and $W_{t,13}^{+}$. However, as the wave age increases in figures 8(e–h) and 8(i–l), $W_{t}^{+}$ is mainly governed by $W_{t,11}^{+}$. This arises due to the strengthening of SRNS and wave-induced velocity gradient by a fast wave (see the supplementary material), thus facilitating wave–turbulence exchange.

Figure 8. Energy transfer between turbulence and waves (wave–turbulence exchange). The total and components ($W_{t}^{+},W_{t,11}^{+},W_{t,13}^{+},W_{t,33}^{+}$) are plotted for (a–d) slow wave (S1 with $ak=0.13$ and $c/u_{\tau }=3.69$), (e–h) intermediate wave (I1 with $ak=0.13$ and $c/u_{\tau }=17.25$) and (i–l) fast wave (F1 with $ak=0.13$ and $c/u_{\tau }=35.80$). The red dashed line is the critical height.

4.3 Correlation between energy transfer and wave regimes

We further correlate the energy transfer with the wave regimes. Firstly, the averaging of $W_{t}^{+}$ leads to $\bar {W}_{t}^{+}=-{{\overline {\langle {{{{u}'_i}}}{{{{u}'_j}}}\rangle {{{\tilde {S}}}_{ij}}}}^{\,+}}=-{{\overline {\widetilde {{{{{u}'_i}}}{{{{u}'_j}}}}{{{\tilde {S}}}_{ij}}}}^{\,+}}$ and $\bar {T}_{t}^{+}=\bar {M}_{t}^{+}+\bar {W}_{t}^{+}$ (detailed information can be seen in the supplementary material), which means the energy transferred into turbulence can be divided into ensemble-averaged and wave-induced parts. Similarly, $\bar {T}_{w}^{+}=\bar {M}_{w}^{+}+\bar {W}_{w}^{+}$ denotes that the energy transferred into wave motion includes the terms of ensemble-averaged production and transportation.

We integrated the profiles of $\bar {M}_{t}^{+}$, $\bar {M}_{w}^{+}$ and $\bar {W}_{t}^{+}$ along the vertical direction to obtain the local-averaged energy transfer. We also defined the total energy as $TE=\bar {M}_{t}^{+}+\bar {M}_{w}^{+}+\bar {W}_{t}^{+}$. The normalized energy components can be expressed as $L_{M}^{t}={\int {\bar {M}_{t}^{+}\,{\rm d}\zeta }}/{\int {TE\,{\rm d}\zeta }}$, $L_{M}^{w}={\int {\bar {M}_{w}^{+}\,{\rm d}\zeta }}/{\int {TE\,{\rm d}\zeta }}$ and $L_{W}^{t}={\int {\bar {W}_{t}^{+}\,{\rm d}\zeta }}/{\int {TE\,{\rm d}\zeta }}$, respectively. Figure 9(a) shows the normalized energy components vary as a function of the wave age. As the wave age increases, the normalized mean TKE production decreases, while the normalized WKE production shows the reverse trend. The normalized wave–turbulence exchange is slightly enhanced with the rising wave age. It is seen that the WKE production contributes negatively to the total energy for the slow wave condition. This energy deficit is balanced through the mean TKE production enhancement. However, an intermediate wave partitions the total energy into nearly equal WKE and mean TKE productions, potentially demonstrating the intermediate wave age as a transition, repartitioning kinetic energy. Then, the WKE contributes more to the total energy than the mean TKE for a fast wave. Notably, the wave–turbulence exchange in the present study is positive and thus can be added into the mean TKE production term and merged as total TKE production. We show $L_{M}^{t}+L_{W}^{t}$ and $L_{M}^{w}$ in figure 9(b). A transitional variation of TKE and WKE productions at intermediate wave age can be clearly observed.

Figure 9. (a) Normalized energy components. (b) Normalized total TKE (TKE produced by wave-coherent and wave-induced motions) and WKE productions.

4.4 Wave growth rate

Following Cao & Shen (2021), the wave-induced pressure can lead to a form drag on the wave surface. This form drag would result in the growth of wave, usually characterized as wave growth rate parameter $\beta$, which is defined as

(4.4)\begin{equation} \beta = \frac{2{{F}_{p}}}{{{(ak)}^{2}}},\quad \text{where}\ {{F}_{p}} = \int_{{0}}^{{{\lambda }_{x}}}{{{{\tilde{p}}}^{ + }}\frac{\partial \eta }{\partial x}\,{{\rm d} x}}. \end{equation}

Here the $+$ denotes normalization by $\rho u_{\tau }^2$. Figure 10 shows the wave growth rate parameter $\beta$ varying with wave age. For comparison, we plot some experimental and numerical results from Plant (1982), Mastenbroek et al. (1996), Cohen (1997) and Cao & Shen (2021). It is seen that the wave growth rate falls into the region overlapping with previous studies, which further verifies the present numerical model. Moreover, a transition of $\beta$ from positive to negative can be found near the intermediate waves.

Figure 10. Wave growth rate varying with the wave age.

5. Concluding remarks

The present paper has investigated the turbulent flow over propagating waves with different wave ages by solving the full Navier–Stokes equations in a moving coordinate frame. After carefully verifying the numerical model, phase-averaging was conducted to extract the ensemble-averaged, wave-induced and turbulent fields and then to emphasize the momentum and energy transfer mechanisms of wind waves.

The careful verification of the numerical model with previous experiments and numerical simulations suggests that simplifying wind waves as turbulent flow over moving boundaries is rational to some extent, even for strongly forced cases (strong coupling of wind waves or slow wave regimes). The dominant mechanism of momentum flux and kinetic energy transfers is related to the wave age. We find that the RSS and TKE productions for slow waves can be equivalent to the two-dimensional shear turbulence with a positive strain rate. As wave age increases, a fast wave, however, enhances the SRNS and windward side's negative RSS, as well as the gradient of both streamwise and vertical velocities. This results in additional productions of RSS and TKE that can be ignored under the slow wave regimes. Nevertheless, these additional productions play a negative contribution, namely suppressing the turbulence. The intermediate wave can be regarded as a transitional condition determining this change.

The energy transfers between wave-coherent, wave-induced and turbulent motions are complex. Specifically, on the leeward side of a slow wave, kinetic energy is transferred from wave-induced to wave-coherent motions and then to turbulence. Meanwhile, wave-coherent motion injects kinetic energy into both wave-induced and turbulent motions on the windward side. For fast waves on the leeward side, turbulence is consumed and transferred into wave-coherent kinetic energy, facilitating wave-induced motion. On the windward side, wave-coherent motion transfers energy into turbulence and suppresses wave-induced motion. The wave–turbulence exchange describes an energy transfer between wave-induced and turbulent motions. Increasing wave age would immensely facilitate this exchange due to the amplification of SRNS and gradient of wave-induced velocity.

The present study reveals in detail how kinetic energy transfers between wave-coherent, wave-induced and turbulent motions, potentially conducive to predicting kinetic energy repartition at different wave regimes. However, there remain some unresolved aspects. Since wind waves, especially for growing conditions, are strong coupling processes, a two-phase air–water interface flow model is more suitable for simulation, which should be constructed in future studies. Furthermore, it should be acknowledged that the flow dynamics of wind–wave interactions is affected by the Reynolds number, which merits further investigation in future work.