Wave–Current Interaction: A 2DH Model for Turbulent Jet and Bottom-Friction Dissipation

Abstract

A correct representation of the non-linear interactions between waves and currents is one of the key points when studying the morphological evolution of nearshore environments, in particular close to river mouths or tidal inlets. Undoubtedly, the numerical modelling of similar phenomena can be very complex and computationally demanding, given the size of the domains. In the present paper, a two-dimensional horizontal (2DH) numerical model is applied to investigate the hydrodynamics of a turbulent jet current interacting with frontal waves, preparatory to the study of morphodynamical processes. The purpose is to reproduce accurately the turbulence of the current flow, which develops in both vertical and horizontal planes, even with the simplifications of depth-averaged velocities. Moreover, the bottom shear stress induces a mechanism of dissipation, which acts both on the jet hydrodynamics and on the wave field. Significant attention is given to this process, which turns out to be crucial in shallow waters. The present model, based on classic shallow-water equations and wave action balance, is applied to a literature test. Comparisons with theoretical and numerical outcomes are shown, the latter obtained with a quasi-3D model.

1. Introduction

The morphological evolution of a river delta or a tidal inlet is the result of the interaction between jet-like current dynamics, tides and wind wave-induced processes. Sediment transport mechanisms, resulting from this interaction, have important effects on the equilibrium configurations of the bottom of channels and of the near inlet region. These processes are governed by the non-linear interaction between currents and gravity waves, one of the main focuses of maritime and coastal hydraulic research.

Theoretical and experimental studies on waves propagating with and against a current have been widely conducted [1,2,3,4], highlighting the effects on wave kinematics (Doppler effect and changes in the wave number and frequency due to shoaling and refraction) and dynamics (changes in wave steepness and wave-action conservation). In turn, the current field is affected by the mutual interaction, through forces generated by the radiation stress tensor of wave fluctuations and increased bed shear stress. In this regard, the bottom friction coefficient of currents can be enhanced significantly by waves [5,6,7,8].

In the nearshore region, characterized by shallower water depth, the interaction between currents and waves with the bottom is one of the main processes that affects the hydrodynamic field.

A river- or tidal-induced current, debouching into the sea, can be regarded as a plane turbulent jet characterized by longitudinal velocity decay and an expansion of its cross section [9]. Several authors have investigated the velocity distribution of a turbulent jet. One of the main contributions is by Abramovich [10], who described the basic principles of a free flow, referring to Prandtl’s hypothesis of free turbulence on the horizontal plane. This theory establishes the form of the similarity function and thus the velocity profile in the jet cross section [11]. Moreover, Öszoy and Ünlüata [12] examined with an analytical approach a turbulent jet, which exits from tidal inlets. They considered several factors: the lateral entrainment, the bathymetric changes and, in particular, the role played by bottom friction. The authors started from widely used depth-averaged equations, thus becoming a reference also for comparison with experimental or numerical results.

Giger et al. [13] and Adami and Milan [14] discussed, from both theoretical and experimental points of view, the effects of bottom friction on the reduction of the longitudinal velocity. They found that this phenomenon had no influence on the crosswise distribution. This outcome is in contrast with the analytical results of Öszoy and Ünlüata [12]. Ismail and Wiegel [15] are the first authors who investigated the behavior of a shallow plane turbulent jet in the presence of opposing surface gravity waves, pointing out the role of the radiation stress components on the flow field.

In nearshore environments, waves have a dominant role that considerably affects the morphodynamics [16,17]. It follows that the correct prediction of the wave field is fundamental to a proper estimate of the morphodynamic changes near the coast.

A key aspect in the spectral transformation of waves, as they propagate from deep to shallow water, is the dissipation of wave energy caused by interaction with the bottom. In particular, the friction coefficient has been observed to vary significantly as a function of both the wave Reynolds number and the relative roughness; the expressions proposed for the calculation of the friction coefficient are numerous and different [18], but precise indications on how to apply them are not always provided. These different formulations can be used in similar contexts, but they can give diverse results, both in the wave field and in the bottom shear-stress values, depending on the wave period. In fact, the wave period may not remain constant in the nearshore propagation, in particular if waves interact with currents. The known Doppler effect can be more or less pronounced according to the angle of incidence between waves and current, for instance in the two cases of an ebbing flow from a river mouth or a tidal inlet, and a longshore current. With the period, the relative roughness also changes and, consequently, the amount of wave-energy dissipation in the presence of shallow depths. The particular aspect of the variability of the wave period and the consequent effects on the bottom friction dissipation in the wave–current interaction, has not been previously dealt with.

It is clear at this point that the bottom strongly conditions the non-linear mechanisms of the wave–current interaction, considerably affecting the current and the wave field, and the bottom shear stress, with important effects on the morphological evolution of littorals or river and tidal inlets. The complexity of the phenomenon necessarily requires dedicated numerical modelling able to investigate hydrodynamic and morphodynamic processes involving estuarine and tidal environments.

A numerical approach needs the development of an appropriate model that should be computationally efficient and able to solve complex hydrodynamic interactions with relative simplicity. In tidal and coastal domains, a widely used methodology consists of coupling different modules [19,20,21,22,23,24]: a hydrodynamic model based on classical shallow water equations with 1D, 2D or a quasi-3D approach and a spectral wave model that reproduces nearshore wave fields, taking into account all the input and dissipation-energy terms.

It has been largely established that two-dimensional horizontal (2DH) models are useful to study the morphological evolution of littorals or nearshore bottoms requiring large domains of application. Moreover, 2DH models can reproduce the observed coastal flow field satisfactorily [22,25].

A quasi-3D model is undoubtedly more complete, but at the same time the improvement in the estimation of vertical eddy viscosity involves a greater computational effort mainly due to the following aspects: the vertical discretization of the domain, which can become demanding in very extensive domains such as coastal or tidal environments, and the calibration process of related parameters. In particular, the last condition is more pronounced when considering the interaction with the vertical profile of wave orbital velocity and radiation stresses, the distribution of which is an on-going debate [26,27]. Furthermore, the logarithmic description of the mean current velocity profile was experimentally found to be valid also in its interaction with waves [28], confirming the validity of a 2DH approach.

Following these considerations and due to the very good hydrodynamic results obtained in literature, Petti et al. [29] presented and applied an efficient 2DH model to a jet current expanding into a shallow basin in the presence of opposing frontal waves. The authors obtained good agreement with both theoretical and numerical results computed by means of a quasi-3D model [30]. The role played by radiation stress forces and the increasing bottom shear stresses in the wave–current interaction was examined [31]. In particular, near the jet inlet, these two mechanisms seem to have the same effect on the flow field, while at greater distances, the contribution to the decay of the longitudinal current velocity due to bottom shear stress is dominant.

In the present paper and starting from these results, the authors investigate the effects induced by bottom shear stresses on the spectral propagation of the waves against the jet current. As indicated above, the bottom dissipation is crucial to the correct prediction of both the wave kinematics and the flow field. In several applications the choice of spectral wave dissipation is not adequately motivated, and sometimes even omitted. For these reasons, the authors intend to verify how the relative roughness affects the wave characteristics, the current velocities, and the bottom shear stress due to the spatial variability of the relative wave period so as to improve the comprehension of the interaction phenomenon, also in the perspective of a morphodynamic application.

Therefore, different approaches used in the computation of the spectral wave dissipation by bottom friction are examined and compared to each other.

In Section 2 the numerical scheme is briefly presented, and in Section 3 the benchmark test is described. Section 4 reports a discussion of the numerical results obtained.

2. Numerical Scheme

In this work, a bidimensional hydrodynamic model has been coupled with the spectral model SWAN (Simulating WAves Nearshore) [32].

The hydrodynamic model solves the classical shallow-water equations:

$$\frac{\partial U}{\partial t}+\nabla \cdot \left[F,G\right]-\nabla \cdot \left[F_t,G_t\right]=S$$

(1)

with:

$$U={\left[h,Uh,Vh\right]}^T;\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}S={\left[0,-gh\frac{\partial z_b}{\partial x}-\frac{{\tau }_{mx}}{\rho }+\frac{F_x}{\rho },-gh\frac{\partial z_b}{\partial y}-\frac{{\tau }_{my}}{\rho }+\frac{F_y}{\rho }\right]}^T;$$

(2)

$$F={\left[Uh,U^{2}h+\frac{g{h}^2}{2},UVh\right]}^T;\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{F}_t={\left[0,2{\upsilon }_{t}h\frac{\partial U}{\partial x},{\upsilon }_{t}h\left(\frac{\partial U}{\partial y}+\frac{\partial V}{\partial x}\right)\right]}^T;$$

(3)

$$G={\left[Vh,UVh,V^{2}h+\frac{g{h}^2}{2}\right]}^T;\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{G}_t={\left[0,{\upsilon }_{t}h\left(\frac{\partial U}{\partial y}+\frac{\partial V}{\partial x}\right),2{\upsilon }_{t}h\frac{\partial V}{\partial y}\right]}^T$$

(4)

(t, x, y) being temporal and horizontal spatial coordinates; U the variable vector; (F, G) and (F_t, G_t) the vectors of advective and turbulent fluxes; and S the source term vector. Moreover, h is the water depth; (U, V) the mean velocities over the depth in x- and y- directions; g the gravity acceleration; $\rho$ the water density; z_b the bottom height; and υ_t the horizontal eddy viscosity coefficient. (F_x, F_y) are the forces due to wave radiation stresses, which depend on the radiation stress tensor S_ij as:

$${\left[F_x,F_y\right]}^T=-{\left[\frac{\partial S_{xx}}{\partial x}+\frac{\partial S_{yx}}{\partial y},\frac{\partial S_{xy}}{\partial x}+\frac{\partial S_{yy}}{\partial y}\right]}^T.$$

(5)

(τ_mx, τ_my) are the components along x and y, respectively, of the mean bottom shear stress τ_m, due to the interaction of current and wave motion. Among all the theories available in literature which combine their effect, in the present paper the formula proposed by Soulsby [31] has been preferred [33]:

$${\tau }_m={\tau }_c\left[1+1.2{\left(\frac{{\tau }_w}{{\tau }_c+{\tau }_w}\right)}^{3.2}\right]$$

(6)

Here, the stresses due to the current τ_c and to the waves τ_w are evaluated as:

$${\tau }_c=\rho \frac{g{n}^2\left[{\left(Uh\right)}^2+{\left(Vh\right)}^2\right]}{h^{7/3}};\hspace{0.17em}\hspace{0.17em}{\tau }_w=\frac{\rho f_w{U}_w{}^2}{2};\hspace{0.17em}\hspace{0.17em}{f}_w=1.39{\left(\frac{A}{k_n/30}\right)}^{-0.52};\hspace{0.17em}\hspace{0.17em}A=\frac{U_{w}T}{2\pi }$$

(7)

with n the Manning coefficient; U_w the maximal bottom orbital velocity; T the relative wave peak period; A the horizontal semi-orbital excursion at the bottom; k_n the equivalent roughness length; and f_w the wave friction factor [31].

The horizontal eddy viscosity coefficient υ_t is evaluated by means of a Smagorinsky approach, as:

$${\upsilon }_t=c_s{}^2{l}_c{}^2\sqrt{2{\left(\frac{\partial U}{\partial x}\right)}^2+2{\left(\frac{\partial V}{\partial y}\right)}^2+{\left(\frac{\partial U}{\partial y}+\frac{\partial V}{\partial x}\right)}^2}$$

(8)

with c_s a numerical coefficient and l_c a characteristic length.

The structure of (8) is similar to υ_t as evaluated following Prandtl’s theory, assuming the product c_s·l_c as the mixing length. In this context, Abramovich [10] found that the mixing length in a free turbulent jet is proportional to the half width of the jet itself b, which varies linearly along the longitudinal axis. In this way, c_s can be considered as the proportional constant and l_c as the half width of the jet. Aiming to adopt a magnitude of b, representative of the whole jet field, its value measured at the inlet b₀ has been chosen. The coefficient c_s has been calibrated in previous tests, based on the experiments carried out by Ismail and Wiegel [15], starting from some values proposed by Abramovich. Finally, c_s has been assumed equal to 0.2 and l_c equal to b₀.

The numerical integration of Equation (1) is carried out by means of a well-balanced second order accurate finite volume method, based on the Harten-Lax- van Leer-Contact (HLLC) Riemann solver associated with a hydrostatic variable reconstruction that assures the scheme is well balanced also in wet and dry conditions [34,35,36].

The wave parameters are evaluated through SWAN [32], an open source finite difference model that solves the wave action density balance equation:

$$\frac{\partial N}{\partial t}+\frac{\partial c_{x}N}{\partial x}+\frac{\partial c_{y}N}{\partial y}+\frac{\partial c_{\vartheta }N}{\partial \vartheta }+\frac{\partial c_{\sigma }N}{\partial \sigma }=\frac{S_{tot}}{\sigma }$$

(9)

where N is the wave action density, defined as E/σ, and it is preferred over wave energy E since it is conserved in the presence of a current field [3]. The propagation velocities in Cartesian and spectral spaces (c_x, c_y, c_σ, c_ϑ), can be written as:

$$\left(c_x,c_y\right)=\frac{dx}{dt}=c_g+U_c=\frac{1}{2}\left[\left(1+\frac{2kh}{\sinh 2kh}\right)\right]\frac{\sigma k}{k^2}+U_c$$

(10)

$$c_{\sigma }=\frac{\partial \sigma }{\partial t}=\frac{\partial \sigma }{\partial h}\left(\frac{\partial h}{\partial t}+U_c\cdot \nabla h\right)-c_{g}k\cdot \frac{\partial U_c}{s}$$

(11)

$$c_{\vartheta }=\frac{d\vartheta }{dt}=-\frac{1}{k}\left(\frac{\partial \sigma }{\partial h}\frac{\partial h}{\partial m}+k\cdot \frac{\partial U_c}{\partial m}\right)$$

(12)

σ being the relative radiant frequency; ϑ the spectral wave direction; k the wave number vector; c_g the wave group velocity vector; U_c the current velocity vector (U, V); and (s, m) a local system of coordinates, respectively s in the direction ϑ and m perpendicular to s. S_tot, the term on the right side of (9), is the source term that takes into account all physical processes including energy generation by wind, dissipation through whitecapping, bottom friction and depth-induced wave breaking, and wave–wave non-linear transfer.

Equation (9) can also represent phenomena like shoaling or refraction, which are related to non-linear interaction and dissipation during wave propagation in shallow water and in the presence of changing water levels and current fields. Among the theories implemented in SWAN to evaluate the dissipation due to bottom friction, in the present paper the approaches by Collins [37] and Madsen et al. [38] are considered, which express the source term due to bottom friction with the same structure and a linear dependence by a bottom-friction coefficient C_f. Collins [37] assumes it as constant and equal to 0.015, while in Madsen et al. [38] C_f is a function of a wave friction factor f_wM related to the equivalent bottom roughness height:

$$C_f=\frac{f_{wM}}{\sqrt{2}}; \{\begin{array}{cc}\frac{1}{4\sqrt{f_{wM}}}+{\log }_{10}\left(\frac{1}{4\sqrt{f_{wM}}}\right)=-0.08+{\log }_{10}\left(\frac{A}{k_n}\right)& \mathrm{if} A/k_n\ge 1.57\\ f_{wM}=0.3& \mathrm{if} A/k_n<1.57\end{array}$$

(13)

SWAN and the hydrodynamic module interact through currents, radiation stress and bottom shear stress. In the present work they have been coupled at discrete time steps (∆t_int). This means that they run separately one after the other, so that wave-radiation stress, wave parameters (U_w and T), and current fields are updated and interchanged with each other.

For each run, the water level and current velocities resulting from the hydrodynamic model are set as input data in SWAN and, on the other hand, the forces due to wave-radiation stresses, the maxima of the orbital motion near the bottom, and the wave relative peak period resulting from SWAN are set as input data in the hydrodynamic model.

Each module works on its own computational grid and the results, to be transferred from one module to the other, are interpolated on an exchange grid.

The whole process is managed by a main program that runs hydrodynamic and SWAN modules alternatively, as suggested by the flow chart depicted in Figure 1, until the simulation ends.

3. Simulation Setup

The test case is that proposed by Petti et al. [29] in comparison with Nardin et al. [30]: a jet issuing in a rectangular basin and interacting with an opposing wave field. As depicted in Figure 2, three different computational grids are built: one mesh for the hydrodynamic model, one for SWAN, and one for the data exchange.

The hydrodynamic mesh is slightly different from that of Petti et al. [29]: it consists of a rectangular channel river inlet entering into a 3000 m large and 2000 m long basin (Figure 2). The elements are rectangular, with a size of 5 m in the x-direction (along the longitudinal axis of the jet) while in the y-direction (along the cross section of the jet) the cell dimension linearly increases from 5 m in a band close to the river centerline up to 50 m in the farthest part of the domain. This guarantees a higher density of grid cells in the central region, where the main characteristics of the flow field develop. Different grid resolutions have been tested, finding differences in the hydrodynamic numerical results in the order of a few percentage points. For this reason, the numerical results can be considered independent of grid size.

The basin has a constant depth of 3 m. A constant discharge of 3.8 m³/sm is assigned to the boundary located at x = −200 m, along a width of 100 m, which coincides with the initial width of the jet. Along the boundary near the jet inlet, a wall boundary condition is assigned, while a mean sea water level coinciding with the still water level is imposed on the remaining sides of the basin.

The resulting velocity downstream of the inlet is 1.4 m/s, corresponding to the velocity indicated in Nardin et al. [30]. A constant Manning coefficient of 0.018 s/m^1/3 is adopted for the hydrodynamic computational mesh, coherent with the value of Chezy coefficient declared by Nardin et al. [30].

The SWAN cartesian computational grid is regular with squared elements of 10 m sides. It has the same water depth but it is larger than the hydrodynamic domain to keep the lateral boundaries sufficiently away from the central area where the wave–current interaction is investigated. Three configurations of waves are imposed at the open sea boundary of the domain, located at x = 2000 m: a “bin wave spectrum”, one for each simulation, is specified, locating energy in one frequency bin closest to a corresponding peak period of respectively 3 s, 5 s and 8 s (Figure 2). The chosen values of the peak period are typical of wind waves generated in fetch-limited and shallow-depth basins, such as the Northern Adriatic Sea and its lagoons. The same significant wave height, equal to 0.5 m, is assigned in all the performed simulations of wave-current interaction. In the remaining sides of the SWAN mesh, outflow boundary conditions were imposed (i.e., waves could leave the area freely).

SWAN is set in non-stationary mode to follow the temporal evolution of the interaction between the waves against the jet current. The physics commands are set up in order to reproduce the steepening induced by currents [1,2,3,4], that can lead to white-capping or, in the limited condition given by the shallow depth, to depth-induced breaking; this last term is modelled in the form proposed by Battjes and Janssen with a constant breaking parameter equal to 0.78 [39].

The wave decay by bottom friction is taken according to the formulations of Collins [37] and Madsen et al. [38] as described in Section 2. Using the Madsen approach, an equivalent roughness-length scale of the bottom of 0.025 m is taken. Petti et al. [29] set the dissipation term for bottom friction conforming to the empirical JONSWAP model for wind-sea conditions [40].

The lowest and the highest discrete frequencies used in the calculation are respectively set equal to 0.05 Hz and 0.595 Hz. The spectral directions cover the full circle with a resolution of 6°.

For better comprehension, the list of the performed simulations with corresponding acronym and characteristics is reported in

Table 1. List of the simulations performed.

Test	Model	Hs (m)	Tp (s)	Spectral Bottom Friction	Wave–Current Interaction
jet alone	hydrodynamic	-	-	-	-
T3C	coupled	0.5	3	Collins	complete 1
T3M	coupled	0.5	3	Madsen	complete
T5	coupled	0.5	5	No friction	complete
T5C_rds	coupled	0.5	5	Collins	radiation stress only
T5C	coupled	0.5	5	Collins	complete
T5M	coupled	0.5	5	Madsen	complete
T5J	coupled	0.5	5	JONSWAP	complete
T8C	coupled	0.5	8	Collins	complete
T8M	coupled	0.5	8	Madsen	complete

¹ Waves and current interact both through radiation stress and the mean bed shear stress τ_m.

.

The last computational grid, named the “exchange grid”, is the smallest and represents the part of the domain in which the hydrodynamic and the spectral model interact. It is a regular grid with squared elements of 10 m sides for a total length of 2000 m in both plane directions. The domain, inside which the jet develops, coincides with that of Petti et al. [29] and Nardin et al. [30].

The interaction between waves and current is performed alternating hydrodynamic and spectral modules every ∆t_int = 300 s until reaching a steady-state condition.

4. Results and Discussion

The hydrodynamics of the jet current alone, without the interaction with a superimposed wave field, is well interpreted by the numerical model [29]. The assumption of mean velocity and, thus, of a logarithmic profile in the vertical section, subtended by the hypothesis of fully developed turbulent flow, is suitable to interpret the behavior of a river or a tidal flow debouching into a large shallow basin.

Following Abramovich [10], two main properties characterize a free turbulent jet current in the main region or otherwise called ZOEF (zone of established flow): the exponential decay of longitudinal velocity and the self-similarity profile of the velocity in cross sections. In shallow depth, bottom friction enhances the velocity decrease, whereas the self-similar profile is maintained. In this regard, Öszoy and Ünlüata [12] proposed the analytical expression for the longitudinal velocity for different values of bottom friction and suggested faster spreading rates of the jet.

In the present test, the hydrodynamic results of the jet alone are compared with the results obtained in the analogous test carried out by Petti et al. [29] and Nardin et al. [30] by means of a quasi-3D numerical model that computes the vertical eddy viscosity with a standard k-ε closure scheme and the horizontal eddy viscosity with a large eddy simulation technique. In Figure 3 the normalized longitudinal velocity along the jet centerline is reported.

The very small differences compared to the previous test [29] could be ascribed to the different boundary conditions used in the present paper. In order to maintain a constant level inside the domain, Petti et al. [29] inserted a crosswise weir on the opposite side of the hydrodynamic mesh compared to the jet inlet, while in the present simulation a fixed water level has been imposed. This last solution is an improvement because it requires a less extensive domain and its applicability is supported by the results obtained. The comparison is also made with the theoretical profile of Öszoy and Ünlüata [12] defined by means of an analytical approach based on depth-averaged shallow water equations, confirming the 2DH model herein used.

Similar considerations can be made about the normalized velocity profile, taken along a cross section and depicted in Figure 4; the two theoretical self-similar profiles proposed by Schlichting and Tollmienn for a free turbulent jet in the main region [10], are shown for comparison with numerical results.

The graph shows that the self-similar profile is maintained along transverse transects but with a greater spread due to bottom friction, according to Öszoy and Ünlüata [12]. This statement justifies the increased width at the base of the section of the numerical curves instead of theoretical ones without friction.

Figure 3 and Figure 4 show good agreement of the present model with numerical quasi-3D and analytical results [29], with the great advantage that it is computationally less expensive than a quasi-3D model.

Petti et al. [29] analyzed the effects of an opposing wave field imposed on the previously described hydrodynamic steady-state condition. The same test is herein repeated, assigning a wave bin spectrum with a peak period of 5 s at the open sea boundary of the domain (Figure 2).

In order to verify the effects of the bottom friction in the wave propagation and in the subsequent interaction with the jet current, diverse tests are performed using different formulations of spectral energy dissipation, as summarized in Table 1.

The most important hydrodynamic evidence compared to the previous simulation is a more pronounced decay of the longitudinal velocity, as depicted in Figure 5, where velocity linear contours are compared in the two cases of the jet flow alone and the jet interacting with waves. For both conditions, the present current field is completely analogous to that reported by Nardin et al. [30].

The flow field develops on the horizontal plane with the prevalent component of the current velocity aligned with the jet axis and a lateral spreading due to the transverse transfer of momentum, which involves formation of eddies on the jet sides. A characteristic feature of a turbulent jet, coherent with both theory and experimental results [10], is the smallness of the transverse velocity components in any section of the jet, compared to the longitudinal ones. In the cases analyzed, the maximum value reached by the y-directed transverse velocity is about 7 cm/s, corresponding to 5% of the velocity at the inlet. For this reason, the waves propagating against the jet current are weakly refracted and remain substantially collinear, at least in the domain portion where the interaction with the current velocities is more important, as pointed out in Figure 6a.

With respect to this, the major changes in the wave field are observed in a narrow region around the jet centerline, where the current vorticity is greater. As can be seen in Figure 6b,c, the main evidence is the steepening of the waves as they approach the jet inlet, where current velocities are higher, due to the increase of wave height and a simultaneous progressive reduction of the relative period as predicted by the Doppler effect. The obtained results are completely analogous to those found by Nardin et al. [30].

Given the underlined peculiarities of both the wave and the flow field, the analysis focuses on the main variables extracted along the jet centerline, starting from the normalized longitudinal velocity of the jet, depicted in Figure 7 and Figure 8. In Figure 9 the normalized longitudinal velocity is taken along a cross section.

As can be appreciated from the graphs, the spectral friction scheme adopted is important and is fundamental for the prediction of the bottom shear stresses, which are also responsible for sediment transport mechanisms. This aspect, in several numerical applications of nearshore wave–current interaction [20,21,22], is not specified.

The major evidence compared to the simulation without wave–current interaction, is the more pronounced decay of longitudinal velocity of the current caused by the increase of the global bottom shear stress due to the wave motion, which leads to the definition of τ_m [31].

In particular, two aspects emerge from the relative comparison of the results: the first is a larger velocity reduction obtained by the present model compared to Petti et al. [29] and Nardin et al. [30] numerical results. Petti et al. [29] have already registered a higher steepening of waves to that calculated by Nardin et al. [30], who do not specify the frictional scheme setup.

However, the present model also gives differences using the same spectral bottom friction approach of JONSWAP. This outcome is attributable to the different equivalent roughness height adopted in the present simulation, suitable to represent a sandy bottom with ripples, consistent with literature references [18,31]. In fact, in nearshore shallow depths, bed shear stresses induce ripples which enhance the bottom dissipation due to form drag. The equivalent roughness height enters directly into the wave bottom shear-stress calculation and consequently into the mean total shear stress that reduces current velocity.

The second important aspect, relevant to this work, is the different longitudinal velocity decay obtained in T5C and T5M simulations. It can be seen in Figure 8 that the discrepancies, due to the bottom-friction effects, increase when moving toward the jet inlet. For low gradients in the velocity of the jet current that occur away from the inlet, the curves substantially coincide. On the other hand, where higher gradients of velocity are registered the two frictional approaches give different results. This happens from the inlet up to x/b₀ = 20. Within this initial region, the Doppler effect of waves against the current is more pronounced, leading to a reduction in the relative wave period as attested also by Nardin et al. [30]. Together with the period, the relative roughness A/k_n also decreases.

Following Kamphuis [41], the bottom-friction coefficient C_f is a function of both the wave Reynolds number, Re = AU_w/ν, where ν is the kinematic viscosity coefficient, and the relative roughness A/k_n. A wave rough turbulent motion near the bed, in which C_f depends only on the relative roughness, can be distinguished from a laminar or transitional one. Hasselmann and Collins [42], using wave spectra measured offshore of Panama City, have indicated a constant value of C_f equal to 0.015. This value has been widely used, also in recent modelling, even though it has a great variability, experimentally verified [18], as found by Madsen et al. [38], who obtained the expression (13) as a function of the relative roughness. Figure 10a,b show that, as the relative roughness decreases, the bottom friction coefficient considerably increases.

Figure 11 shows values of the Reynolds number and the relative roughness obtained in the present simulations performed. These values are compatible with a wave rough turbulent flow at the bottom, confirming the validity of the bottom-friction approaches applied.

However, in Madsen formulation, the spectral energy dissipation is increased by the progressive steepening of the waves against the jet current and consequently the enhancement of the bottom-friction coefficient. This explains the higher longitudinal velocity in the T5M simulation compared to T5C (Figure 8).

In the mechanisms of interaction between the waves against current, a gradient of radiation stress tensor also arises from the steepening of the waves and, thus, additional hydrodynamic forces enter into the motion equations.

Petti et al. [29] have already highlighted that in the first region of the spreading jet, up to x/b₀ = 10, the contribution to the decay of the velocity due to radiation stress and bottom shear stress is similar, due to the presence of major gradients of radiation stress toward the inlet. For this purpose, the simulation T5C has also been performed taking into account only the radiation stress contribution and neglecting the wave-bottom shear stress in the computation of τ_m. The results are very close to that of Nardin et al. [30], as can be seen in Figure 7.

To check if a different wave period can further influence the jet hydrodynamics, simulations with, respectively, a peak period of 3 s and 8 s have also been carried out. In analogy with the above, the results in term of normalized longitudinal current velocity are reported in Figure 12 and Figure 13.

The comparison between the simulations with the different assigned periods does not highlight significant differences in the jet flow field; the curves relating to the same bottom friction model are practically overlapping. However, the same cannot be said if the wave bed shear stress (Figure 14a) and the maximum derived from the wave–current interaction (Figure 14b) are investigated.

According to Soulsby [31], the maximum bed shear stress is expressed as follows:

$${\tau }_{\max }={\left[{\left({\tau }_m+{\tau }_w\cos \varphi \right)}^2+{\left({\tau }_w\sin \varphi \right)}^2\right]}^{1/2}$$

(14)

$\varphi$ being the angle between the current velocity and the direction of wave travel. This quantity is fundamental in morphodynamics because it is responsible for the lifting of sediments from the bottom when it exceeds the threshold bed shear stress for the onset of motion.

If the wave bed shear stresses are considered, a slight dependence on the period is noticed. Increasing the periods, the wave bed shear stresses are higher and reach the maximum value closer to the jet inlet. Moreover, comparing the results with the same period, great differences are found related to the bottom friction scheme, consistent with the previous results.

The maximum bottom shear stress, deriving from the wave–current interaction, is significantly affected by greater values of the wave period. The configurations performed with a peak period of 3 s do not show particular differences either among themselves or with the bottom shear stress of the jet alone. The contribution of wave shear stress to the maximum bottom shear stress is small for low periods at the assigned depth. Instead, τ_w related to 5 s and 8 s peak period, increase τ_max considerably and differently if the bottom friction coefficient C_f is left to vary with the relative bottom roughness [38].

This outcome can be very important in modelling hydrodynamic nearshore processes, especially in shallow water depth environments, where the interaction with bottom plays a crucial role, both on currents and waves, as can be seen from the numerical results described above.

For this purpose, in a context, for instance, close to a tidal inlet or an estuarine environment, spectral energy dissipation presents great variability related to the peak period, which can change significantly in the presence of strong gradients of the velocities due to the Doppler effect. This mechanism depends on the angle between current and waves: the Doppler effect is not as evident in the case of quasi-orthogonal incidence, as occurs for the longshore currents and the waves propagating nearshore, but it becomes important if current and waves are quite collinear, as for river or tidal jets. For these reasons, it is difficult to consider that the frictional coefficient can be regarded as constant over the whole domain.

The effects of this dissipation are particularly expressed on the bottom shear stresses, responsible for sediment resuspension and motion and the morphological evolution of the bottom itself.

5. Conclusions

In the present work, a simple 2DH numerical model has been applied to a literature test [29,30] reproducing a river jet current issuing into a large shallow basin and interacting with frontal waves. The model couples a hydrodynamic numerical model with the spectral model SWAN and it is suitable for reproducing the interaction between waves and currents typical of estuarine and coastal areas. Moreover, the known computational advantages of 2D modelling allow it to be applied to large domains.

The outcomes have been compared to both analytical and numerical results, the latter being obtained with a model that adopts a more complex representation of the vertical and horizontal eddy viscosity. The comparison shows that the present model provided a very adequate representation of the hydrodynamics of a turbulent jet confined to the bottom.

In detail, the role played by bottom friction has been investigated and it has been recognized as the dominant process in the decay of the jet longitudinal velocity and its greater spreading. Bottom friction dissipation becomes even more important in the correct reconstruction of the nearshore wave field, above all in shallow-depth contexts.

The simulations of the jet current interaction with opposing travel waves carried out have shown that, in the spectral mechanism of dissipation energy, the proven dependence of the bottom friction coefficient with the relative roughness at the bed can change the wave field significantly. Consequently, the interaction with current can also affect hydrodynamics.

To this end, two different theories implemented in SWAN to evaluate the dissipation due to bottom friction have been considered: the Collins approach [37], which assumes the bottom friction coefficient as constant, and the Madsen scheme [38], which takes into account the variability of the same coefficient in contexts of wave rough turbulent motion at the bed.

The effects of these two approaches are quite sensitive to the period of waves, which can change in the domain due to shoaling and Doppler mechanisms generated by currents. In particular, bottom friction has a dominant role in the calculation of the maximum shear stress, which is responsible for bed-erosion processes. The stress is significantly increased by greater values of the wave period and of the current velocity. Above all, it is conditioned by the assumed spectral dissipation scheme.

From a general point of view, both of the theories analyzed are suitable for calculating energy dissipations that wind waves, propagating on shallow domains, undergo due to interaction with the bed. However, near a river delta or a tidal inlet, where currents can be strong and quite collinear compared to the mean wave direction, relative roughness can vary, consequently affecting the friction coefficient and the related bottom shear stress. The Madsen formulation [38] takes this variability into account, as opposed to the Collins scheme [37], and for this reason it results as preferable for interpreting the interaction phenomenon.