Numerical and Experimental Investigation of Wave Overtopping of Barriers

: We present a study of wave overtopping of barriers. The phenomenon of the wave overtopping over emerged structures is reproduced both numerically and experimentally. The numerical simulations are carried out by a numerical scheme for three-dimensional free-surface ﬂows, which is based on the solution of the Navier–Stokes equations in a novel integral form on a time-dependent coordinate system. In the adopted numerical scheme, a novel wet–dry technique, based on the exact solution of the Riemann problem over the dry bed, is proposed. The experimental tests are carried out by adopting a nonintrusive and continuous-in-space image-analysis technique, which is able to properly identify the free surface even in very shallow waters or breaking waves. A comparison between numerical and experimental results, for several wave and water-depth conditions, is shown.


Introduction
The prediction of the wave overtopping process over barriers is fundamental to properly design seawalls, breakwaters, sea dikes and, more generally, structures that have the aim to protect inland areas from sea waves. For this reason, a large number of studies (e.g., [1][2][3][4][5]) has been carried out in the past years with the aim of giving estimation to the fundamental parameters that govern the aforementioned phenomenon. These studies led to the development of modern overtopping prediction formulae [6]. This approach provides a method of evaluation of several essential parameters (overtopping discharge, maximum run-up, wave transmission), being a robust and simple tool for the prediction of the key aspects of wave overtopping.
An alternative strategy to study the wave overtopping phenomenon consists in the use of numerical models. In order to have a reliable prediction of the wave overtopping over a structure, a numerical model has to consistently simulate different phenomena: wave transformation from deep to shallow water, wave breaking over variable bathymetry, wave run-up on the structure, wave transmission, three-dimensional effects. The first numerical models presented in literature for the overtopping simulation were based on solution of the depth-averaged Nonlinear Shallow Water Equations (NSWE) (e.g., [7][8][9]). Using shock-capturing methods, these NSWE-type models possess the capability to naturally simulate wave breaking [8]. Unfortunately, in NSWE models, frequency dispersion effects are neglected, so wave propagation cannot be properly reproduced. To overcome this A set of experiments used to validate the numerical results, has been carried out in the wave flume of the Hydraulics Laboratory of the Department of Civil and Environmental Engineering and Architecture, University of Cagliari, Italy. Field [27] and laboratory measurements are widely used to validate numerical models [28,29]. Many of these measurements are carried out by means of resistive probes as experimental instruments, by means of which water level can be measured. Using resistive probes, the prediction of the water level in very shallow water zones and in the wave-breaking zone, is a challenging issue. In this context, Lorke et al. [3] used thin gauges and micropropellers to predict the water level and flow velocity at the crest of the structure. In order to overcome the aforementioned issue, we developed a nonintrusive and continuous-in-space technique, based on Image Analysis, developed and tested by Ferrari et al. [30]. This technique is able to properly identify the wave free surface even in prohibitive situations for traditional resistive probes, such as very shallow waters and/or breaking waves. Furthermore, a validation of the proposed model with a random wave overtopping test, is presented; in fact, random wave tests are of great practical interest, because they can represent complex real phenomena that can be found in nature.
In this work, we present a numerical and experimental study for the wave overtopping over barriers. The main element of novelty in the adopted numerical model, consists in the implementation of wet-dry technique by means of which the celerity of the wet-dry front is correctly evaluated even over flat or negative slope bottoms. The paper is structured as follows: in Section 2, we describe the numerical model; in Section 3, we describe the setup for the experimental tests; in Section 4, we present the results of several validation tests for the numerical model; in Section 5, we present the conclusion of the study.

Governing Equations
Let us consider a moving control volume ∆V(t), with surface ∆A(t), which moves with a velocity where ρ is the density of the fluid, → f = −grad(Gx 3 ) is the external body force vector per unit mass (where G is the constant of gravity), T is the stress tensor,n is the unit outward normal vector and the symbol represents the outer product. For an incompressible flow, T can be defined with the constitutive relation T = −PI + 2µS (where P is the total pressure, S is the strain rate tensor, I is the identity tensor and µ is the dynamic viscosity). In the context of a free-surface flow, we define η as the free-surface elevation, with respect to an horizontal reference plane. The total pressure can be split into an hydrostatic part ρG(η − x 3 ) and a dynamic part q. Therefore, Equation (1) can be rearranged as: where ν = µ/ρ is the kinematic viscosity. The integral from of the continuity equation over the moving Assuming that the flow is incompressible, Equation (3) becomes: Let (ξ 1 , ξ 2 , ξ 3 ) be a system of curvilinear coordinates. The transformation from Cartesian coordinates (x 1 , x 2 , x 3 ) to the generalized curvilinear coordinates (ξ 1 , ξ 2 , ξ 3 ) is The metric tensor is c lm = → c (l) · → c (m) , with (l, m = 1, 2, 3). The Jacobian of the transformation is given by Let us consider ∆V(t) as a volume element defined by surface elements bounded by curves lying on the coordinate lines. We define the volume element in the physical space as: and the volume element in the transformed space as: It is possible to see that the volume element defined in Equation (7) is time dependent, while the one defined in Equation (8) is not. Similarly, we define the surface element in the physical space as ∆A(t), which is time dependent, and the surface element in the transformed space as ∆A * , which is constant over time.
In a generalized curvilinear coordinate system Equations (2) and (4) become, respectively: and In Equations (9) and (10) and hereinafter, ∆A * α+ and ∆A * α− indicate the contour surfaces of the volume element ∆V * on which ξ α is constant and which are located at the larger and the smaller value of ξ α , respectively, and α, β, γ = 1, 2, 3 are cyclic.
In order to simulate the fully dispersive wave processes, Equations (9) and (10) can be transformed in the following way. Let x 3 = 0 at the horizontal reference plane defined by the still free surface. Let where h is the still water depth. Let the bottom elevation with respect to the horizontal reference plane be z(x 1 , x 2 ) = −h(x 1 , x 2 ). Our goal is to accurately represent the bottom and surface geometry and correctly assign the pressure and kinematic conditions at the bottom and at the free surface. A particular transformation from Cartesian to curvilinear coordinates, in which coordinates vary in time in order to follow the free-surface movements, is Under the transformation (Equation (11)), the components of vector This coordinate transformation basically maps the time-varying coordinates of the physical domain into a fixed coordinate system (ξ 1 , ξ 2 , ξ 3 ) where ξ 3 spans from 0 to 1. In addition, the Jacobian of the transformation becomes √ c = H We define the cell-averaged value (in the transformed space), respectively of the conservative variable H → u and of the primitive variable H (recalling that H does not depend on ξ 3 ): where ∆A * ξ 1 ξ 2 = ∆ξ 1 ∆ξ 2 is the horizontal surface element in the transformed space. By recalling that ∆V * is not time dependent, by dividing by ∆V * and by recalling Equations (13) and (14), Equation (9), expressed in the time-dependent coordinate system defined in Equation (11), becomes: By recalling that ∆A * is not time dependent, by dividing by ∆A * ξ 1 ξ 2 , and by recalling Equations (13) and (15), Equation (10), expressed in the time-dependent coordinate system defined in Equation (11), becomes: (17) in which ξ α+ and ξ α− indicate the contour lines of the surface element ∆A * on which ξ α is constant and which are located at the larger and the smaller value of ξ α , respectively. Equation (17) represents the governing equation that predicts the free-surface motion. Equations (16) and (17) represent the expression of the three-dimensional motion equations as a function of the H → u and H variables in the coordinate system (ξ 1 , ξ 2 , ξ 3 ). The numerical integration of the abovementioned equations allows the fully dispersive wave-propagation simulation. In order to Water 2020, 12, 451 6 of 17 take into account the effects of turbulence, we introduce a turbulent kinematic viscosity, estimated by means of the Smagorinsky subgrid model.

Numerical Procedure
Equations (16) and (17) are solved by means of the numerical model proposed by [26]. The numerical scheme consists in a combined finite-volume and finite-difference scheme with a Godunov-type method, and it is briefly described as follows.
A staggered grid is used, in which velocity is evaluated at cell centers, while pressure is evaluated at horizontal cell faces. Cell-averaged values H → u and H are known at the center of the calculation cells and are the known variables. t (n) is the time level of the known variables while t (n+1) is the time level of the unknown variables. A three-stage, third-order nonlinear strong-stability-preserving (SSP) Runge-Kutta scheme is used in the advancing of the solution.
In order to obtain a divergence-free velocity field at each time level, we apply a pressure correction formulation. Knowing Hu l (n) , we evaluate Hu l (n+1) using the three-stage iteration procedure described below. Let Hu l (0) = Hu l (n) At each stage p (with p = 1, 2, 3) an auxiliary field Hu l (p) * is obtained from Equation (16) by using values from the previous stage where D(H, u l , τ) indicates the right-hand side of Equation (16), in which the term related to the dynamic pressure gradient is neglected. See [31] for the values of coefficients Ω pq , ϕ pq and d q . In order to numerically solve Equation ( l * from Equation (19) and by calculating H (p−1) * from Equation (17), does not satisfy the continuity equation. To obtain a nonhydrostatic divergence-free velocity fields, we correct the pressure field and the auxiliary velocity field, at each stage p, by introducing a scalar potential Ψ, to be evaluated by the well-known Poisson-like pressure equation, Equation (20) is approximated by means of a second-order cell-centered finite-difference scheme, so that it can be reduced to an algebraic linear system with the form: where A is the coefficient matrix (with 15 nonzero diagonal coefficients), Ψ is the scalar potential vector and b is the vector of constant terms. The algebraic linear system (Equation (21)) is solved by means of an iterative procedure based on a four-color, zebra-line Gauss-Seidel alternate method and a multigrid V-cycle accelerator [32]. The corrector irrotational velocity field is evaluated as follows: In order to have a divergence-free velocity field at each stage, we correct the velocity field as follows: Let us indicate with L(H, u l , τ) the right-hand side of Equation (17). The advancing of the depth to the p stage is obtained by: Hu l (n+1) is given by In order to simulate the run-up and the backwash dynamics of the wet and dry front on a structure, a novel wet-dry numerical technique is adopted in the scheme. This technique allows the scheme to simulate precisely the advance of the wet-dry front on the structure.

Wet-Dry Technique
In the adopted numerical model, a novel wet-dry technique is implemented, in order to allow the model to simulate waves propagating over flat or negative slope bottoms. Previous research [15,26] adopted a simple technique proposed by [24] for the treatment of the wet-dry front movement. In order to introduce the aforementioned technique, let us define a column as a set of calculation cells stacked on top of each other, having the same x 1 and x 2 coordinates. In [15,24,26], a dry column changes its status to wet if its free-surface elevation is less than the free-surface elevation of an adjacent wet column. The aforementioned treatment is proven to be effective in the case in which the wet-dry front propagates over a positive sloping bottom (see Figure 1a), i.e., when at the wet-dry front the bottom elevation of the dry column is greater than the one of the adjacent wet column. Otherwise, if the bottom elevation of the dry column is equal to or less than the one of the adjacent wet column (flat or negative sloping bottom, Figure 1b), the technique proposed by [24] leads to an instantaneous wetting of the subsequent dry columns with bottom elevation equal to or less than the one of the first wet column. This would lead to a thin-film treatment, in which the dry state is treated as a wet state with a thin layer of water; as stated by [33], the thin-film technique fails to correctly evaluate the celerity of the wet-dry front.
In the context of the study of the wave overtopping over barriers, the bottom can have a zero or negative slope. In order to avoid the aforementioned drawback, the following wet-dry procedure has been implemented in the numerical scheme. For the sake of simplicity, the technique is described for the x 1 direction. Let (i, k) be the two indexes which define the center of the calculation cells in the x 1 and x 3 directions, respectively. Hereinafter, the superscript (n) indicates the time level τ n of the known variables, while the superscript (n + 1) indicates the time level τ n+1 = τ n + ∆τ, of the unknown variables. Let us define the column I i as the set of calculation cells stacked on top of each other, having the same i index. The water depth in a wet column is always greater than a minimum water depth, H i > H min . In a dry column H i = H min . Let I i be a dry column and I i−1 a wet column. The criterion by means of which I i changes status from dry to wet is described as follows.
Reconstructions defined at the cells that form the wet column I i−1 lead to point-value variables i−1/2 , the exact solution of the wet-dry Riemann problem over a dry bed is used, adopting piecewise constant initial data, H ID and u ID : (29) in which x 1i−1/2 represents the position in the x 1 direction of the interface between the columns I i−1 and I i . From [33], we have:S Is to be noted that the standard solution of the Riemann problem over wet bed, with nonzero water depth in the initial data, cannot be used, since it would lead to a wrong celerity evaluation. See [33] for the details about the wet-dry Riemann problem solution.
Knowing the wet-dry front celerity from Equation (30), the distance between x 1i−1/2 and the wet-dry front is: The column I i becomes wet if the wet-dry front reaches the successive column interface, so that where ∆x 1 is the calculation cell dimension in the x 1 direction. The above described procedure is repeated in the x 2 direction. It can be noted that, by using the proposed wet-dry technique, the celerity of the wet-dry wave front is used to evaluate the time in which the subsequent dry cell becomes wet. In this way, in case of a flat or negative sloping bottom, a dry column does not instantaneously change status to wet.
Water 2020, 11, x FOR PEER REVIEW 8 of 17 Is to be noted that the standard solution of the Riemann problem over wet bed, with nonzero water depth in the initial data, cannot be used, since it would lead to a wrong celerity evaluation. See [33] for the details about the wet-dry Riemann problem solution.
Knowing the wet-dry front celerity from Equation (30), the distance between ⁄ and the wet-dry front is: The column becomes wet if the wet-dry front reaches the successive column interface, so that where ∆ is the calculation cell dimension in the direction. The above described procedure is repeated in the direction. It can be noted that, by using the proposed wet-dry technique, the celerity of the wet-dry wave front is used to evaluate the time in which the subsequent dry cell becomes wet. In this way, in case of a flat or negative sloping bottom, a dry column does not instantaneously change status to wet.

Experimental Setup
The laboratory tests to validate the numerical model were carried out in a 21.00 m long, 0.30 m wide and 0.50 m high flume, with glass walls, a piston-type wave-maker and an absorbing beach designed to minimize the reflections. The wavemaker is controlled by an in-house developed software, so that monochromatic regular wave trains are produced. The software controls the wave period and amplitude. See Figure 2 for the details of the experimental setup. More details on the wavemaker and flume can be found in [30] and [34].
The model breakwater employed in the numerical simulations was reproduced by means of a black-painted trapezoidal obstacle of Perspex. Three hydrostatic water depths were employed (0.29,

Experimental Setup
The laboratory tests to validate the numerical model were carried out in a 21.00 m long, 0.30 m wide and 0.50 m high flume, with glass walls, a piston-type wave-maker and an absorbing beach designed to minimize the reflections. The wavemaker is controlled by an in-house developed software, Water 2020, 12, 451 9 of 17 so that monochromatic regular wave trains are produced. The software controls the wave period and amplitude. See Figure 2 for the details of the experimental setup. More details on the wavemaker and flume can be found in [30,34].

Numerical Simulation Details
All the numerical simulations were run on a 12-CPU workstation using the Message Passing Interface (MPI) procedure. In all the numerical simulations, 75%-80% of the computational time was spent for the solution of the Poisson-like equation in the velocity correction step of the numerical procedure. For this reason, a multigrid V-cycle accelerator was used in the numerical model (see Section 2.2), in order to reduce the simulation time.

Three-Dimensional Validation Test
In order to validate the three-dimensional nonhydrostatic numerical model presented, we reproduced the experimental test performed by [37] on a fixed submerged barrier with periodically spaced rip channels.
According to [37], the basin is characterized by a length equal to 17.20 m and width equal to 18.20 m. The beach has an initial steep slope (1:5) followed by a milder slope (1:30). There are three   The model breakwater employed in the numerical simulations was reproduced by means of a black-painted trapezoidal obstacle of Perspex. Three hydrostatic water depths were employed (0.29, 0.26 and 0.23 m). The water was seeded with a fluorescent dye, the investigation area was lighted with a light sheet and images of the experiments were recorded by a digital video camera placed orthogonal to the investigation area; consequently, the recorded images showed a bright area (water) and a dark area (background and breakwater). A set of runs with the same combination of wave periods (0.95 and 1.55 s), wave heights (3.50, 2.50 and 1.00 cm) and hydrostatic water levels of the numerical tests was performed.
In order to not modify the flow and to overtake the issues related to traditional probes with very shallow water or breaking waves, an image-analysis technique was developed to identify the free surface. Image-analysis techniques [34,35] have many advantages when compared to traditional probes, as they are nonintrusive (i.e., the flow is not altered by the probes) and quasicontinuous in space (i.e., a measure can be obtained on every pixel on the recorded images). A review can be found in [36].

Numerical Simulation Details
All the numerical simulations were run on a 12-CPU workstation using the Message Passing Interface (MPI) procedure. In all the numerical simulations, 75-80% of the computational time was spent for the solution of the Poisson-like equation in the velocity correction step of the numerical procedure. For this reason, a multigrid V-cycle accelerator was used in the numerical model (see Section 2.2), in order to reduce the simulation time.

Three-Dimensional Validation Test
In order to validate the three-dimensional nonhydrostatic numerical model presented, we reproduced the experimental test performed by [37] on a fixed submerged barrier with periodically spaced rip channels.
According to [37], the basin is characterized by a length equal to 17.20 m and width equal to 18.20 m. The beach has an initial steep slope (1:5) followed by a milder slope (1:30). There are three submerged breakwaters: the two lateral ones are 3.66 m long, and the central one is 7.32 m long. The distance between the submerged breakwaters is 1.82 m. The seaward edges of the bar sections are located at approximately x = 11.10 m with the bar crest at x = 12.00 m and their shoreward edges at x = 12.30 m. In Figure 3, a plain view and a longitudinal section of the basin for the experimental test performed by [37] are shown. edges at = 12.30 m. In Figure 3, a plain view and a longitudinal section of the basin for the experimental test performed by [37] are shown.
Regular monochromatic waves are generated at the offshore boundary of the domain. The test parameters are deep water wave height = 5.12 cm, wave period = 1 s, average water depth at the bar crest ℎ = 4.73 cm and cross-shore location of the still water line = 1490 cm. The computational grid resolution is ∆ = 0.025 m, ∆ = 0.05 m and six layers are used in the vertical direction. The time step is ∆ = 0.0025 s. The simulation was carried on for a simulated time of 365 s, and ran for four days. After 65 s, the hydrodynamic field is supposed to be quasistationary, so the time-averaging of the hydrodynamic quantities is carried out in the interval 65 s-365 s.   Figure 4 shows a comparison of the cross-shore currents at two different longshore sections: the first at x = 12.20 m (onshore side of the bar) and the second at x = 13.00 m (between the bar and the shoreline). From Figure 4, it can be seen that at the onshore side of the bar (x = 12.20 m), the agreement between the numerical results and the measurements obtained by [37] is good; from Figure 5, it can be seen that in the onshore zone (x = 13.00 m), the numerical results are in quite good accordance with the experimental ones, even if the numerical model slightly overestimates the offshore currents near the gap between the barriers.
Water 2020, 11, x FOR PEER REVIEW 10 of 17 Figure 4 shows a comparison of the cross-shore currents at two different longshore sections: the first at = 12.20 m (onshore side of the bar) and the second at = 13.00 m (between the bar and the shoreline). From Figure 4, it can be seen that at the onshore side of the bar ( = 12.20 m), the agreement between the numerical results and the measurements obtained by [37] is good; from Figure 5, it can be seen that in the onshore zone ( = 13.00 m), the numerical results are in quite good accordance with the experimental ones, even if the numerical model slightly overestimates the offshore currents near the gap between the barriers.
In order to show the ability of the model to reproduce three-dimensional velocity fields, the numerical model is validated against the laboratory measurements of [38], obtained in the same test configuration described above. Figure 5 shows the vertical distribution of time-averaged cross-shore normalized velocity profiles. Normalization is done by dividing the velocity by the celerity = √ ℎ. The results are obtained at = 9.00 m, = 13.60 m (2.00 m offshore of the bar) and at = 11.75 m, = 13.60 m (inside the gap between the bars). Numerical results are compared against the experimental measurements obtained by [38]. From Figure 5, it can be seen that the agreement between the numerical results and experiments is very good. The model is able to reproduce both the offshore current near the free surface in the offshore zone (Figure 5a) and the strong offshore current located at mid-depth in the zone inside the gap between the bars (Figure 5b).    . Cross-shore currents at x = 12.20 m (a) and at x = 13.00 m (b). Circles: experimental results by [37]. Solid line: proposed numerical model.

Wave Overtopping Test
To predict the behavior of waves interacting with the barrier in different configuratio experimental tests were carried out (see Table 1). In test T1, mean water depth is lower than barrier height, resulting in an emerged barrie T2, mean water depth is equal to barrier height, resulting in a partially emerged barrier. In mean water depth is higher than barrier height, resulting in a submerged barrier. For results phase-averaged free-surface elevation is computed. Numerical tests are carried out w proposed scheme to reproduce the experiments. For the numerical tests, a fixed time step 0.0004 s, a horizontal spatial step of ∆ = 0.01 m and six layers in the vertical direction are u simulations were carried on for a simulated time of 600 s, and ran for approximately 24 hou 40 , the hydrodynamic field is supposed to be quasistationary, so the phase averaging of surface elevation is carried out in the interval 40 s-600 s. Figure 6 shows the comparison between numerical results and experimental measurem test T1. It must be noted that the wave is able to pass over the barrier only partially. Figure 6 the wave at the end of the run-up phase; the run-up celerity of the wave in the numeric slightly lower than the measured one. Figure 6b shows the wave propagation over the cre barrier; numerical and experimental results have a good agreement, suggesting that the p numerical model is well able to predict front propagation over flat bottoms, thanks to the p wet-dry technique. Figure 6c shows that even the simulated wave run-down is in agreem the measured one. Figure 6d shows the wave trough located at the offshore side of the barrie Figure 6e shows the beginning of the run-up phase; Figure 6d,e displays a good accordance numerical and experimental predictions. Note that in Figure 6c,d, at the onshore side of the some disturbances of the experimental measurements are present. In order to show the ability of the model to reproduce three-dimensional velocity fields, the numerical model is validated against the laboratory measurements of [38], obtained in the same test configuration described above. Figure 5 shows the vertical distribution of time-averaged cross-shore normalized velocity profiles. Normalization is done by dividing the velocity by the celerity c = √

Test T1
Gh. The results are obtained at x = 9.00 m, y = 13.60 m (2.00 m offshore of the bar) and at x = 11.75 m, y = 13.60 m (inside the gap between the bars). Numerical results are compared against the experimental measurements obtained by [38].
From Figure 5, it can be seen that the agreement between the numerical results and experiments is very good. The model is able to reproduce both the offshore current near the free surface in the offshore zone (Figure 5a) and the strong offshore current located at mid-depth in the zone inside the gap between the bars (Figure 5b).

Wave Overtopping Test
To predict the behavior of waves interacting with the barrier in different configurations, three experimental tests were carried out (see Table 1). In test T1, mean water depth is lower than barrier height, resulting in an emerged barrier. In test T2, mean water depth is equal to barrier height, resulting in a partially emerged barrier. In test T3, mean water depth is higher than barrier height, resulting in a submerged barrier. For results analysis, phase-averaged free-surface elevation is computed. Numerical tests are carried out with the proposed scheme to reproduce the experiments. For the numerical tests, a fixed time step of ∆t = 0.0004 s, a horizontal spatial step of ∆x = 0.01 m and six layers in the vertical direction are used. The simulations were carried on for a simulated time of 600 s, and ran for approximately 24 h. After 40 s, the hydrodynamic field is supposed to be quasistationary, so the phase averaging of the free-surface elevation is carried out in the interval 40 s-600 s. Figure 6 shows the comparison between numerical results and experimental measurements for test T1. It must be noted that the wave is able to pass over the barrier only partially. Figure 6a shows the wave at the end of the run-up phase; the run-up celerity of the wave in the numerical test is slightly lower than the measured one. Figure 6b shows the wave propagation over the crest of the barrier; numerical and experimental results have a good agreement, suggesting that the proposed numerical model is well able to predict front propagation over flat bottoms, thanks to the proposed wet-dry technique. Figure 6c shows that even the simulated wave run-down is in agreement with the measured one. Figure 6d shows the wave trough located at the offshore side of the barrier, while Figure 6e shows the beginning of the run-up phase; Figure 6d,e displays a good accordance between numerical and experimental predictions. Note that in Figure 6c,d, at the onshore side of the barrier, some disturbances of the experimental measurements are present.  Figure 7 shows the comparison between numerical results and experimental measurements for test T2. In this test, the wave reaches the crest of the barrier with a very shallow water depth. Figure  7a shows the wave run-up over the offshore side of the barrier; as in the test T1, is to be noted that the run-up is well predicted and that numerical and experimental results have a good agreement. Figure 7b,c shows that the proposed wet-dry technique allows the numerical model to reproduce well the wave propagating over the crest of the structure. In Figure 7d, the wave is shown propagating over the inshore side of the barrier crest; numerical and experimental results are in a good agreement, even if is to be noted that experimental results show that the offshore side of the crest is dry. In Figure 7e, the run-down over the offshore side of the barrier is shown, in which numerical and experimental results both well predict the phenomenon.  Figure 7 shows the comparison between numerical results and experimental measurements for test T2. In this test, the wave reaches the crest of the barrier with a very shallow water depth. Figure 7a shows the wave run-up over the offshore side of the barrier; as in the test T1, is to be noted that the run-up is well predicted and that numerical and experimental results have a good agreement. Figure 7b,c shows that the proposed wet-dry technique allows the numerical model to reproduce well the wave propagating over the crest of the structure. In Figure 7d, the wave is shown propagating over the inshore side of the barrier crest; numerical and experimental results are in a good agreement, even if is to be noted that experimental results show that the offshore side of the crest is dry. In Figure 7e, the run-down over the offshore side of the barrier is shown, in which numerical and experimental results both well predict the phenomenon. Figure 7b,c shows that the proposed wet-dry technique allows the numerical model to reproduce well the wave propagating over the crest of the structure. In Figure 7d, the wave is shown propagating over the inshore side of the barrier crest; numerical and experimental results are in a good agreement, even if is to be noted that experimental results show that the offshore side of the crest is dry. In Figure 7e, the run-down over the offshore side of the barrier is shown, in which numerical and experimental results both well predict the phenomenon.  Figure 8 shows the comparison between numerical results and experimental measurements for test T3. In this test, the barrier is completely submerged. Figure 8a,b shows the wave crest propagating over the offshore side of the barrier, with good agreement between numerical and experimental results. In Figure 8c,d, the wave-front propagation over the structure is well reproduced; in the experimental prediction the free surface is slightly higher than in the numerical one. Figure 8e shows the wave crest passing the barrier; is to be noted that in the numerical simulation the wave front has a slightly higher celerity than the predicted experimental celerity.

Test T3
Water 2020, 11, x FOR PEER REVIEW 13 of 17 Figure 8 shows the comparison between numerical results and experimental measurements for test T3. In this test, the barrier is completely submerged. Figure 8a,b shows the wave crest propagating over the offshore side of the barrier, with good agreement between numerical and experimental results. In Figure 8c,d, the wave-front propagation over the structure is well reproduced; in the experimental prediction the free surface is slightly higher than in the numerical one. Figure 8e shows the wave crest passing the barrier; is to be noted that in the numerical simulation the wave front has a slightly higher celerity than the predicted experimental celerity.

Random Wave Overtopping Test
In this subsection, we show the results of several validation tests (described in Table 2), initially proposed by [39], for the proposed numerical model, in the context of wave overtopping by random waves. This is of great practical interest, because, in contrast to regular waves, random waves better represent the complexity of the phenomena found in nature.

Random Wave Overtopping Test
In this subsection, we show the results of several validation tests (described in Table 2), initially proposed by [39], for the proposed numerical model, in the context of wave overtopping by random waves. This is of great practical interest, because, in contrast to regular waves, random waves better represent the complexity of the phenomena found in nature. The channel is 80 m long, with a flat bottom in the first 40 m of the channel, and a seawall next to the flat bottom, with a slope angle of 1 : 4. The freeboard is R C = 1 m in the first set of three tests and R C = 1.5 m in the second set of three tests. The still water depth is H SWL = 8 m. Random waves are generated by means of the Jonswap spectrum, with a spectral enhancement factor of γ = 3.3. According to [39], the input random wave parameters are shown in Table 2. For the numerical simulations, a fixed time step ∆t = 0.0002 s, a horizontal spatial step of ∆x = 0.2 m and five layers in the vertical direction are used. The simulations ran for approximately 16 h.
The results obtained by the proposed model, in terms of mean overtopping rate, are compared against the numerical results obtained previously [39,40], and against the results obtained with the wave overtopping formula of [2] (see Figure 9). According to [40], who reproduced the experiment with a Smoothed-Particle Hydrodynamics (SPH) numerical model, the wave overtopping rate becomes stable at approximately t = 50 s. Thus, the mean overtopping rate is calculated in the time interval 50-180 s. The comparison shows that the wave overtopping rate is generally well predicted by the numerical models. The proposed model generally gives a slightly lower prediction than the model proposed by [39]. In the comparison against the prediction formula of [2], the proposed model slightly underestimates the wave overtopping rate in the tests with small wave overtopping, while giving a better prediction in the other tests. In Figure 10, the free-surface elevation and velocity field are shown for test R2 at different instants. The overall shape transformation features of the wave propagating over a slope are shown to be qualitatively well captured.
Water 2020, 11, x FOR PEER REVIEW 14 of 17 Figure 9. Wave overtopping rate. Red: calculated by the Van der Meer and Janssen laboratory formula [2]; green: computed by Soliman [39]; blue: computed by Shao et al. [40]; yellow: computed with the proposed model. Figure 9. Wave overtopping rate. Red: calculated by the Van der Meer and Janssen laboratory formula [2]; green: computed by Soliman [39]; blue: computed by Shao et al. [40]; yellow: computed with the proposed model. Figure 9. Wave overtopping rate. Red: calculated by the Van der Meer and Janssen laboratory formula [2]; green: computed by Soliman [39]; blue: computed by Shao et al. [40]; yellow: computed with the proposed model. The results obtained by the proposed model, in terms of mean overtopping rate, are compared against the numerical results obtained previously [39,40], and against the results obtained with the wave overtopping formula of [2] (see Figure 9). According to [40], who reproduced the experiment with a Smoothed-Particle Hydrodynamics (SPH) numerical model, the wave overtopping rate becomes stable at approximately = 50 s. Thus, the mean overtopping rate is calculated in the time interval 50-180 s. The comparison shows that the wave overtopping rate is generally well predicted by the numerical models. The proposed model generally gives a slightly lower prediction than the

Conclusions
In this work, we presented a numerical and experimental study about the waves overtopping over structures. The numerical model for the simulation of nonhydrostatic free-surface flows proposed by [26] was used and improved, in order to predict the wave overtopping phenomenon. In particular, novel wet-dry technique was implemented in the model, in order to properly calculate the advancement of the wave front over a dry bottom.
The numerical model was validated against several numerical and experimental test cases. Several laboratory tests were carried out, reproducing regular waves overtopping structure, by means of a nonintrusive and continuous-in-space image-analysis technique. Furthermore, in order to evaluate the capacity of the model to simulate phenomena that can be found in nature, the ability to predict the wave and current development over complex structure geometries and the random wave overtopping of structures, were tested. From the numerical tests, it can be assumed that the numerical model is able to reproduce the wave overtopping and wave-front propagation over emerged barriers, for different wave, water-depth and bottom-slope conditions, thanks to the implementation of the proposed wet-dry technique.