Modelling of Wave Energy Dissipation over Porous Mounds using Boussinesq-type Equations

The predictive skills of a Boussinesq-type model which is developed to simulate wave transformation over porous beds are further investigated with a new set of wave flume data. The governing equations of this model include a combined mass conservation equation, momentum conservation equation for the water layer and a momentum conservation equation for the porous medium. Momentum conservation equation for the porous layer incorporates an empirical DarcyForchheimer type term for laminar and turbulent frictional losses and an inertial term. Weakly nonlinear and highly nonlinear wave conditions are generated in the tests to assess the validity of the model parameter values obtained from an optimization technique. The model simulations show excellent agreement with experimental data for weakly nonlinear non-breaking waves and for breaking waves, when the values of the model parameters are determined by the optimization technique.


Introduction
Various methods are being used to mitigate coastal erosion due to both natural and manmade causes around the world. Some resort to soft measures such as beach nourishment, berm management, etc., and some tend to rely on hard measures such as offshore breakwaters, seawalls, revetments etc. Many of these hard measures however involve structures, hence a solid understanding of their interaction with waves is vital.
Coastal structures are mostly designed to be porous in order to induce extra wave energy dissipation by turbulence in the interstices. Therefore, the accurate prediction of wave energy dissipation through porous damping is crucial in producing a better hydraulic design for porous coastal structures. The recent developments in wave theories and numerical models have substantially helped coastal researchers to study the evolution of waves and currents in the vicinity of solid structures but not to the same extent when it comes to porous ones. This is primarily due to the uncertainty in determining a few empirical coefficients associated with porous media.
Flaten and Rigg [1] were the first to propose a Boussinesq-type model for the transformation of water waves over porous bottoms in shallow water with an asymptotic approach. Cruz et al. [2] derived a set of 2DH Boussinesq-type equations for a porous bed of arbitrary thickness and tested their applicability on wave transformation over porous submerged breakwater with an opening for non-breaking waves. They incorporated frequency dispersion into the flow in porous media following Madsen et al. [3]. Hsiao et al. [4] presented a fully nonlinear Boussinesq-type model following the same approach of Nwogu [5] to improve the dispersion properties of the flow inside the porous media. Chen [6] improved Hsiao et al. [4] model by eliminating weak vertical vorticity assumption and retaining higher-order porous damping terms in the momentum conservation equation. The fundamental properties for wave evolution (i.e. wave celerity and porous damping rate) of Chen [6] model equations are investigated in Cruz and Chen [7] for constant depth and horizontal porous interface. Similar research directions with such two-phase approach (i.e., a separate clear water layer and a porous layer) can be found in Johnson et al. [8] and Mettallinos et al. [9]. In all these models, there exist a mass conservation equation , a  momentum conservation equation for clear  water  flow  and  another  momentum  conservation equation for water flow inside  porous media.  The commercial software such as MIKE21 BW,  BOUSS-2D allow waves to penetrate into  porous media by keeping only single  momentum  conservation  equation  for  combined effects, hence cannot be used to  model wave propagation over porous media. On the other hand, the two-phase models will not allow the free water surface to exist in the porous media.
In the present study, the transformation of weakly and highly nonlinear waves over porous mounds is simulated using a model developed based on a truncated form of Chen [6] Boussinesq-type equations. This work is inspired based on the potential of such a model in simulating the wave evolution over submerged breakwaters. This paper primarily discusses the validity of model parameter values (which determine the dispersion and wave damping characteristics of the model), which are obtained using an optimization technique for different wave environments. The paper also presents a new set of data which would help the researchers working along similar directions to validate their numerical models.

2.1
The Governing Equations A model developed based on a truncated form of Chen [6] Boussinesq-type equations is used to simulate the wave energy dissipation over porous beds in the present study. A short description of the governing equations used in the model is presented below.
In one horizontal dimension (1DH), the truncated form of Chen [6] equations (i.e. by dropping higher-order nonlinear and porous damping terms) read as: Continuity equation: Momentum equation for clear water: Momentum equation for porous layer: where, is the porosity of the porous material, is the free water surface elevation, ℎ is the clear water depth, ℎ is the porous layer depth, ℎ =(ℎ + ℎ ) is the total water depth, is the water particle velocity at reference water depth , and is the water particle velocity at reference water depth inside porous layer along the wave propagation direction. The factors and , are determined by matching dispersion characteristics of wave motion in both clear water and porous layer (not discussed in this paper).
is the porous resistance force evaluated at = as follows: where and are respectively, the linear and nonlinear porous resistance coefficients, and is the added mass coefficient. The porous resistance coefficients are estimated from the following relationships following Sollitt & Cross [10]: where is the kinematic viscosity of water (∽ 10 −6 m 2 /s). The intrinsic permeability, , added mass coefficient, , and nonlinear drag coefficient, , are determined by empirical formulae proposed by van Gent [11]. where, 50 is the median diameter of the porous material, , and are empirical coefficients that, in principle, are determined by experiments. However, van Gent [11] concluded that, although the coefficients and depend on the parameters like shape, aspect ratio or the orientation of the porous material, the values 1000 and 1.1 can be used, for them respectively, if the characteristic length scale 50 is used to define the size of porous material. In addition, he suggested using a constant value for empirical coefficient, of 0.34, for simplicity neglecting the complex dependency of coefficient, , on the flow field.

2.2
Energy Dissipation due to Bottom Friction and Wave Breaking Two additional terms are introduced into momentum conservation equation to simulate the dissipation of wave energy due to bottom friction and wave breaking, which are represented by and .
where, is the wave friction factor, and is the near-bottom wave orbital velocity along the x-axis. Following Kennedy et al. [12], a simple eddy viscosity type of formulation is adopted to simulate energy dissipation due to wave breaking.
where, a turbulence closure model is used to determine the eddy viscosity coefficient, . Wave breaking induced energy dissipation is observed in only two cases, hence not discussed in detail in this present paper. Offshore Open Boundary and Moving Shoreline Boundary As there is significant wave reflection from porous submerged mounds (Figure 1), the waves which propagate back to wave incident boundary are absorbed to prevent a build-up of wave energy inside the computational domain using line boundary method proposed by Ishii et al. [13]. Kennedy et al. [12] slot technique is employed to simulate moving shoreline.

2.4
Numerical Scheme The truncated form of Chen [6] Boussinesq-type equations are solved in the time domain using a third-order Adams-Bashforth predictor step and a fourth-order Adams-Moulton corrector step. The computations are carried out primarily for free surface elevation, , reference velocity for clear water, , and reference velocity for porous medium, .

Experimental Set-up
The wave damping characteristics of the model investigated using a new set of experimental data, which include 10 and 02 cases associated with non-breaking and breaking waves, respectively. The set-up includes trapezoidal submerged gravel mounds constructed at the horizontal bed section of a wave flume which is 35.0 m long, 0.6 m wide and 1.0 m deep at the University of Tokyo, Japan. The wave paddle is capable of absorbing any reflected wave preventing building-up of wave energy inside the flume. The porous mounds were constructed with gravel of mean diameter (d50) equal to 20 mm, which had a porosity of 0.44.
The tests were conducted with various incident wave heights, wave periods, offshore water depths, mound heights and widths to represent a wide range of experimental conditions. All cases were tested with regular waves.
are Table 1

-Wave Parameters and Physical Model Dimensions Associated with all Cases
Wave parameters  Figure 1 illustrates the definition of variables associated with the physical models used, and  For waves propagating over an impermeable bed, Nwogu [5] showed that if a proper is used, the range of applicability of the Boussinesq-type equations can be extended to intermediate water depths. The same idea was employed by Hsiao et al. [4] to select a set of and , so that the optimal linear dispersion properties of the model equations can be achieved over a wide range of depths. As Chen [6] Boussinesq-type equations (or truncated forms) are a direct extension of Hsiao et al. [4] equations, they too contain the parameters  and  , which define the vertical locations and of the representative velocities, and , in the water layer and the porous material layer, respectively, as follows:

Model Characteristics
The linear dispersion relation for the truncated form of Chen [6] Boussinesq-type equations is obtained by linearizing Eqs. (1), (10) and (13) as follows: + ℎ + ℎ + 1 ℎ 3 + 2 ℎ 3 = 0 …(22) + + 1 ℎ 2 + 2 ℎ 2 = 0 …(23) where,  is the angular frequency and k is a complex number, where the real part of k denotes the wave number and the imaginary part represents the spatial porous damping rate. In order to solve Eq. (37), a rearrangement of the terms is necessary as follows.   Table 2. The , parameter values are incremented by 0.001 in the optimization process for each case, and the best combination is determined from one million combinations of and  . In cases 01-07, the depth of the porous mound is fixed at 100 mm, and the total water depth is varied from 202 mm to 333 mm, obtaining 03 different relative porous layer thicknesses. In cases 08-12, the depth of porous mound and the total water depth are fixed at 200 mm and 255 mm, respectively, obtaining a constant relative porous layer thickness.

4.2
Model Performance After determining the optimal , values for each case, the model is run to obtain the temporal variation of free water surface elevation over the porous mounds. In cases 01-07, wave gauges were set at 15 locations, and the wave heights were obtained from the records. In cases 08, 10, and 12, wave records were taken at 14 locations and in cases 09 and 11 they were taken at 27 locations (since cases 09 and 11 correspond to breaking wave conditions). Figures 2 and 3 depict the simulated and measured wave height distributions over the porous mounds associated with cases 01-07 and cases 08-12, respectively. The blue lines indicate the wave heights simulated by the model when , parameters are determined from optimization technique (this is called run 1 for convenience hereafter). Since, , parameter values are determined by linearizing the equations in the optimization technique, it is expected that the profiles will be either overestimated or underestimated when the waves become nonlinear. In addition, the representative relative porous layer thickness, , for the computation of , values is set to be constant (taken as the maximum relative porous layer thickness for each case), which is not precise either. Hence, another set of simulations is performed to obtain the values of , parameters considering the best agreement with measured wave height data (this demonstrates the validity of the model) for comparison. In this case, the overall best agreement between measured and simulated values at all gauges is quantified by calculating Willmott index [15] given in Eq. (40) for wave height.
In Eq. (40), denotes the measured values, denotes the computed values, ̅ denotes the average measured values, and is the number of gauges. The Willmott index takes the value of 1.0 for a perfect agreement. The , parameter values which correspond to the overall best agreement between measured and simulated wave heights for each case are listed in Table 2 as well. The red lines in Figures 2 and  3 depict the wave height distributions over the porous mounds corresponding to these simulations (this is called run 2 for convenience hereafter).

Figure 2 -Wave Height Distribution over the Porous Mound Associated with Cases 1-7
When cases 01-07 are considered, the wave heights obtained by run 1 and run 2 agree fairly well with experimental data and indicate only minor deviations from each other (run 2 anyway is the best fit). The Willmott index is found to be in between 0.72 -0.93 and 0.74 -0.95 for run 1 and run 2, respectively (Figure 4). The highest difference in Willmott index between run 1 and run 2 is found to be for case 03 amongst this set of tests. Case 03 is associated with the highest wave nonlinearity parameter amongst these seven cases; hence comparatively more vulnerable to an error in determining , parameters through the optimization technique. The wave nonlinearity parameter, ( /ℎ ), and frequency dispersion parameter, (ℎ/ ), for each case are computed based on the wave height, water depth and wavelength at the centre of the porous mound and are presented in Figure 5. All seven cases are associated with comparatively higher frequency dispersion parameters but are safely within the dispersion properties of the model. Overall, these tests verify that the model is capable of performing well with reasonable accuracy for weakly nonlinear waves, when , parameters are obtained using the optimization technique.

Figure 3 -Wave Height Distribution over the Porous Mound Associated with Cases 8-12
When cases 08-12 are considered, the wave heights obtained by run 1 demonstrate major deviations from the measured data, particularly for case 08 and case 10, as it can be seen from Figure 3. The Willmott indices are found to be in between 0.49 -0.95 and 0.92 -0.97 for run 1 and run 2, respectively ( Figure 4) for this set of tests. The wave nonlinearity is relatively high in these cases; hence, the error due to the linearization in optimization technique is expected to be significant. The frequency dispersion parameter is however found to be less than 0.08 for all five cases as it can be seen from Figure 5. The case 12 which is associated with a regular wave train of period 1.60 s too demonstrates deviations for the wave heights on the porous mound. Therefore, the linearization made to match the dispersion properties is found to be significantly affecting the predictive skills of the numerical model when the waves are highly nonlinear and nonbreaking.

Figure 4 -Willmott Indices Associated with all Cases
The cases 09 and 11 correspond to breaking wave conditions, and the predictive skills of the model are found to be excellent for even run 1. This is rather contradictory considering the high nonlinearity of the waves observed in these two cases. A number of simulation runs conducted with various , combinations demonstrated that the overall agreement is less sensitive to the values of , when the waves are breaking over the porous mound. This is due to the significance of energy dissipation due to wave breaking compared to that of porous damping. This clearly suggests the potential of model's extension to simulate wave evolution over submerged structures nearshore. The Effect of Porosity A separate set of simulations was carried out on the 12 cases taking the porous mound to be impervious to highlight the significance of porous damping in wave evolution. Figure 6 depicts the simulated wave height distributions considering submerged mounds to be impervious (setting = −0.531 following Nwogu [5]) and porous (run 2) associated with cases 04, 07, 09 and 12.
These cases represent two weakly nonlinear wave conditions, one breaking wave condition and one highly nonlinear, non-breaking wave condition. For all cases, a reduction can be observed in wave height distribution due to porous damping with a significant change for case 12 (height of the porous mound is 100 mm for cases 04 and 07 and it is 200 mm for cases 09 and 12).

4.4
Free Water Surface Elevation The temporal variation of the free water surface is compared here (qualitatively) over the porous mounds since the wave height only provides a phase-averaged quantity associated with wave propagation. The case 07 (with weakly nonlinear non-breaking waves), the case 09 (with breaking waves) and the case 12 (with highly nonlinear non-breaking waves) are selected for comparison. Figures 7, 8 and 9 depict the time series of measured and simulated (run 2) free water surface elevations at 15 gauge positions for cases 07, 14 gauge positions for case 09 and case 12 over the porous mounds. The model demonstrates reasonably good predictive skills, capturing the disintegration of higher harmonics in the lee of the porous mound (primarily for cases 09 and 12 where the waves are highly nonlinear over the porous mounds). However, it is noted that the same and values are used for the computation of wave evolution for a given case, although the relative porous layer thickness varies from zero to a finite value within the section of the porous mound. This can lead to a mismatch in both porous damping and frequency dispersion characteristics of the waves within the domain of the porous mound. Furthermore, waves which are approaching and leaving the porous mound propagate only on the water phase. The parameter is invalid for this region and, according to Nwogu [5], the value of the parameter should be equal to -0.531. This can also add an error to the frequency dispersion of waves outside the porous mounds.

Conclusions
The predictive skills of a Boussinesq-type model for simulating wave transformation over porous beds are evaluated with a particular emphasis on the influence of model parameters and on the wave height distribution. The model is found to re-produce the wave damping over porous beds with a high level of accuracy if the two parameters and are set appropriately.
The optimization technique, which is used to determine the values of and parameters for a given relative porous layer thickness, is found to be valid only when the waves are weakly nonlinear. Model simulations demonstrate that the non-breaking waves with high nonlinearity are damping at a different rate when compared with observations. However, highly nonlinear breaking waves are found to be less sensitive to the values and parameters due to the significance of wave breaking induced energy dissipation compared to that of porous damping, which allows the model to be applicable to nearshore. Still, care must be taken when setting values for and parameters offshoreward of the porous mounds, as frequency dispersion will be highly important when the water phase only is present. Therefore, it is highlighted here the necessity of determining the value of and parameters for simulations with an emphasis on wave nonlinearity and changing relative porous layer thickness.