A numerical model of wave overtopping on seadikes

Coastal Engineering 57 (2010) 757–772

Contents lists available at ScienceDirect

Coastal Engineering j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / c o a s t a l e n g

A numerical model of wave overtopping on seadikes Thieu Quang Tuan a,⁎, Hocine Oumeraci b a

Faculty of Marine and Coastal Engineering, Water Resources University (WRU), 175 Tay Son, Dong Da District, Hanoi, Viet Nam Department of Hydromechanics and Coastal Engineering, Leichtweiss Institute for Hydraulic Engineering and Water Resources (LWI), Technical University of Braunschweig, Beethovenstrasse 51a, 38106 Braunschweig Germany

b

a r t i c l e

i n f o

Article history: Received 16 January 2009 Received in revised form 9 April 2010 Accepted 13 April 2010 Available online 15 May 2010 Keywords: Nonlinear shallow water equations Finite volume method Seadikes Surface roller Wave overtopping Overtopping events Combined wave and surge overtopping

a b s t r a c t This paper describes the development of a numerical model for wave overtopping on seadikes. The model is based on the flux-conservative form of the nonlinear shallow water equations (NLSW) solved with a high order total variation diminishing (TVD), Roe-type scheme. The goal is to reliably predict the hydrodynamics of wave overtopping on the dike crest and along the inner slope, necessary for the breach modelling of seadikes. Besides the mean overtopping rate, the capability of simulating individual overtopping events is also required. It is shown theoretically that the effect of wave breaking through the drastic motion of surface rollers in the surfzone is not sufficiently described by the conventional nonlinear shallow water equations, neglecting wave setup from the mean water level and thus markedly reducing the model predictive capacity for wave overtopping. This is significantly improved by including an additional source term associated with the roller energy dissipation in the depth-averaged momentum equation. The developed model has been validated against four existing laboratory datasets of wave overtopping on dikes. The first two sets are to validate the roller term performance in improving the model prediction of wave overtopping of breaking waves. The last two sets are to test the model performance under more complex but realistic hydraulic and slope geometric conditions. The results confirm the merit of the supplemented roller term and also demonstrate that the model is robust and reliable for the prediction of wave overtopping on seadikes. © 2010 Elsevier B.V. All rights reserved.

1. Introduction Inundation and flooding as the result of dike breaching during storms can bring about catastrophic consequences to low-lying coastal areas. There exist dozens of causes of dike breaching, amongst which wave overtopping is the major and becomes increasingly threatening in the present context of global climate change. The predictive capacity of wave overtopping on coastal dikes is therefore essential for risk and safety assessment of hinterland areas. For the dike design at present, the mean rate of wave overtopping is usually sufficient, for which empirical formulae (see e.g. EurOtop, 2007) are commonly used. However, besides the mean rate, the temporal behaviour of wave overtopping is also of interest in many coastal practices such as for drainage design, stability of the inner dike slope, or dike breach modelling. To this end, numerical models appear to be more advantageous and have attracted many researchers in the past decades.

⁎ Corresponding author. Tel.: + 84 4 35634415; fax: +84 4 35635832. E-mail addresses: [email protected], [email protected] (T.Q. Tuan), [email protected] (H. Oumeraci). 0378-3839/$ – see front matter © 2010 Elsevier B.V. All rights reserved. doi:10.1016/j.coastaleng.2010.04.007

Numerical models based on Reynolds-Averaged Navier–Stokes equations (RANS) and on Nonlinear Shallow Water (NLSW) equations are commonly used for the computation of wave overtopping. Use of either of these two model families, hereinafter called RANS models and NLSW models, respectively, should be justified depending on particular purposes of applications. RANS models are much more sophisticated and thus computationally more expensive (they generally require hours in computation time on a normal PC to simulate seconds of flow in real time), but in return flow fields can be described in great detail and for almost any complex bed configurations. Theoretically, full presentation of hydrodynamic processes would be ideal for an advanced processbased breach model, allowing for better descriptions of many complex breaching processes. Practically, RANS models, however, at present are most applicable to problems for which the bed topography is relatively fixed and also the considered time scale is relatively short (up to a few hours only; see e.g. Losada et al., 2008; Reeve et al., 2008). On the other hand, bed changes in a storm-induced breach are rapid and the time scale is longer (on a storm scale from many hours to half a day), which make the application of RANS models at present for breaching of seadikes very challenging and not yet feasible. In contrast to RANS models, NLSW models are much less sophisticated and thus more computationally efficient since only depth-averaged parameters of the flow are computed. Recently, the

758

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

use of NLSW models for studies on wave overtopping on dikes has been very successful. Models of this kind have been applied to reliably predict the overtopping volume even on dikes with very complex geometry and steep seaward slopes (see e.g. Hu et al., 2000; Shiach et al., 2004). Surprisingly, the underlying assumption of hydrostatic pressure distribution of the NLSW equations still seems to work reasonably well in those cases. The present research aims at modelling the hydrodynamics of wave overtopping on dikes, which is needed for the dike breach modelling in the next stage. In essence, it requires the model capacity to reliably predict both the mean overtopping rate and the characteristics of overtopping events. Moreover, the model should be computationally efficient but sufficiently reliable. With this aim, there are a number of reasons that support the use of a NLSW model. First, regarding the hydrodynamics of wave overtopping, the seaward slope of seadikes considered herein is sufficiently gentle so wave overtopping is expected not to be violent like in the case of steep or vertical walls. In addition, the overtopping flow on the dike crest and along the inner slope has been found to reasonably resemble a boundary layer flow without strong turbulence (see e.g. Stansby and Feng, 2004; Schüttrumpf and Oumeraci, 2005). Therefore, the use of the NLSW equations is certainly appropriate for describing wave overtopping on seadikes. Second, at present, our understanding of morphodynamics is generally still lagging far behind that of hydrodynamics. Hence, very detailed flow parameters do not necessarily help improving the model prediction of the breach morphology. For an efficient breach model at present a NLSW-based model might therefore be the better choice. However, at the cost of high efficiency, NLSW models contain several intrinsic drawbacks, reducing the model performance considerably. In this paper, we shall deal with one of these which is the problem of under-predicting the mean water level in the surfzone due to insufficient description of the roller energy dissipation in breaking waves. This paper is organized as follows. In Section 2, the basic model formulation is addressed, including the derivation of a new source term related to the roller energy dissipation for incorporating in the momentum equation. Section 3 covers a discussion of how to incorporate the roller effect in the model. The model validation against four available datasets is presented in Section 4. Finally, conclusions are drawn in Section 5. 2. Numerical model

The nonlinear shallow water (St. Venant) equations (NLSW) in the flux-conservative form read: → ∂U + ∂t

∂x

→
→ = S x; U

ð1Þ

→→ → → → where conserved vectors U, F (x, U), and source term vector S(x, U ) arising from bed slope and friction are defined as follows: → U ðxÞ =

A

Sf =

1 Q jQ j 1 Q jQ j = fw 2 2 gA R C 2 A2 R d

IP = ∫ ðd−zÞBðzÞdz Sbx = −

dZb dx

where C is the Chezy coefficient, fw (= 2g / C2) is the representative (Fanning) friction coefficient, d is the flow depth, R is the hydraulic radius, and B(z) is channel width at an elevation z above the bed, Zb is the bed level. It is noted that for most coastal models B and R are absent from Eqs. (1)–(3). However, it is preferable herein to retain the channel width B for the future use in the modelling of dike breach development (i.e. B is also a variable). To cope with numerical difficulties arising from bed level variations and possible flow discontinuities, a robust shock-capturing numerical scheme should be used. Since exact Riemann solvers are less efficient, approximate Riemann solvers are more preferable. For the present model, a Roe-type Riemann solver in conjunction with the finite volume method (FVM) is adopted (approximate Riemann solver of Roe (1981); see also Toro, 1997). Discretization of Eq. (1) using FVM yields: U

n+1

 Δt
n n =U − −F F + ΔtS Δ x i +1 = 2 i−1 = 2

k 1 1 2 k ˆk ðF + Fi1 Þ− ∑ α j ˆλ jRˆ 2 i 2 k = 1 i1 = 2 i1 = 2 i1 = 2

ð5Þ

where Δ(•)i + 1/2 = (•)i + 1 − (•)i; α̂ki + 1/2, λ̂ki + 1/2 (=û ± ĉ), and Rk̂i ± 1/2 (= (1λ̂ki + 1/2)T ) are the Roe's coefficients (wave strengths), the eigenvalues and the eigenvectors of the so-called Roe's approximate Jacobian, respectively, and û and ĉ are the Roe-averaged velocity and celerity. The 2nd order TVD numerical fluxes are modified from those in Eq. (5) with an additional flux limiter term as follows (see also Dodd, 1998; Hubbard and Dodd, 2002): h

i k 1 1 2 k k ˆk ˆk ðF + Fi1 Þ− ∑ α 1−L ri1 = 2 1−νi1 = 2 jˆλi1= 2 jR i1= 2 2 i 2 k = 1 i1 = 2

ð6Þ

Q

0 1 Q
 → → B C F x; U = @ 2 A Q + gIP A ! 0 →
→
 S x; U = gA Sbx −Sf

ð4Þ

in which n is a known computational level, Δt and Δx are time and space steps, respectively, Fi − 1/2 and Fi + 1/2 are numerical fluxes at the interfaces of a computational cell based on the approximate Riemann solver of Roe (see Toro, 1997). For the sake of consistency some basic formulations are repeated hereafter, reference is made to Toro (1997) for complete details. We present here both 1st order and 2nd order total variation diminishing (TVD) schemes. The 1st order Roe's numerical fluxes read:

Fi1 = 2 =

!

ð3Þ

0

Fi1 = 2 =

2.1. Basic flow formulation

→
→ ∂ F x; U

respectively, IP is a hydrostatic pressure force term acting on the wetted area.

ð2Þ

in which g is the gravitational gravity, A is the wetted cross-sectional area, Q is the flow discharge, Sbx and Sf are bed slope and friction slope,

where νik± 1/2 = λ̂ki ± 1/2Δt / Δ x are local wave Courant numbers and rki ± 1/2 are ratios of the wave strengths. L(r) is a chosen flux limiter function (minmod function was adopted herein), when L(r) = 0 the scheme reduces to the first order one (Eq. (5)). To be consistent with the TVD fluxes, the source term vector is recast in the following form (see Hubbard and Dodd, 2002): Si1 = 2 =

h

i 1 2
k k ˆk ˆ R ∑ 1  sk 1−L ri1 = 2 1−νi1 = 2 β k i1 = 2 2 k=1

ð7Þ

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

where sk = sign(λ̂ki ± 1/2) and β̂ k (k = 1, 2) is a discretized source term coefficient. It is shown in Garcia-Navarro and Vazquez-Cendon (2000) that the Roe scheme is well-balanced if the source term is also upwinded, whereby Eq. (4) now becomes: U

n+1

n

=U −

  Δt 1 Δt
1 ðψL Þi−1 = 2 + ðψR Þi+1 = 2 Fi+1 = 2 −Fi−1 = 2 + Δx Δx 2 2 ð8Þ

where (ψL)i − 1/2 and (ψR)i + 1/2 are source term fluxes contributed from left and right of the nodal point i, respectively. ˆ ψL = β ˆ ψR = β

! s1 fˆ1 −s2 fˆ
2
 ; ˆ 2 1 + s2 fˆ2 ˆ 1 1 + s1 fˆ1 − λ λ ! ˆ −s fˆ s2 f 2 1
1
 ˆ 1 1−s fˆ − λ ˆ 2 1−s fˆ λ 1 1

ð9Þ

2 2

where f ̂k = 1 − L(rik± 1/2)(1 − νik± 1/2). The source term coefficient βk̂ can be determined based on the concept of well-balance scheme when the water level is at rest (see Garcia-Navarro and Vazquez-Cendon, 2000):
 ˆ 1 = −β ˆ2 = β ˆ = gΔx Aˆ Sˆ bx − Sˆ f β ˆ 2c

ð10Þ

in which Â, Ŝbx, and Ŝf are the Roe averages at a computational cell specified according to a point-wise manner (see e.g. Garcia-Navarro and Vazquez-Cendon, 2000). The upwind discretization of the source term vector as shown in Eqs. (8)–(10) is valid for both 1st order and 2nd order TVD schemes presented herein. It is worth mentioning that extra terms in the → → source term vector S (x, U ) in Eq. (2) might be needed in order to account for additional effects in some particular situations, e.g. the effect of wave breaking through the drastic motion of surface rollers in the surfzone as presented in Section 2.3. In such a case, coefficient β̂ in Eq. (10) has to be modified accordingly. Eq. (8) is solved in time to describe the flow conditions, which requires specification of the boundary conditions on both sides of the dike. For the present situation, two boundary conditions, viz. two time series of the water level, one at the upstream inflow boundary and one at the downstream outflow boundary are sufficient (see Fig. 1). Importantly, the inflow boundary should carefully be treated so that it is transparent, avoiding numerical contamination from any reflected waves. Therefore, the present model adopts the characteristics-based approach of Burguete and Garcia-Navarro (2001). In addition, the surface variation of any reflected wave train, needed in the total water level prescribed at the inflow boundary, is calculated from the

759

seaward-advancing characteristics (see also Kobayashi and Wurjanto, 1989). The present model is a wet-bed model, i.e. any dry cells if exist are flooded with an infinitesimal water layer dF of the order of 10− 6 m, while other flow variables of these cells are set to zero. Besides, a cell is regarded as dry if the cell flow depth drops below dF. This procedure allows for the use of the same scheme for determining the interface numerical fluxes. To deal with the problem of wet/dry fronts, it is adopted here the approach by Dodd (1998) whereby the Riemann problem at a wet/dry front within the computational domain is solved using an exact solution from one rarefaction fan (see also Toro, 2001). 2.2. Wave overtopping of breaking waves Most recent numerical studies on wave overtopping have focused on the prediction of average discharge only (see e.g. Dodd, 1998; Hu et al., 2000; Hubbard and Dodd, 2002). However, for many engineering applications as aforementioned, in addition to the average discharge (averaged over the entire storm duration), the characteristics of individual wave overtopping events are also of importance. This study is mainly concerned with wave overtopping on seadikes. A common feature is that seadikes are mostly situated on shallow foreshores and thus intensive wave breaking is expected (see TAW, 2002). However, inspection of the performance of the NLSW modelling of wave overtopping as briefly addressed below reveals that models based on the conventional NLSW equations tend to considerably underestimate the average overtopping rate as well as the number of overtopping events for breaking waves. Kobayashi and Raichle (1994) were amongst the first to address the numerical modelling of overtopping events. The model, called RBREAK, was essentially based on the previous work of Kobayashi and Wurjanto (1992), in which a system of conservation equations of mass, momentum and energy was derived from the two dimensional Reynolds equations. Wave breaking was assumed to occur numerically in the flow field. The used numerical scheme is a finite difference one, which is different from the FVM being used in the present research. Through the model validation against 12 small-scale test runs of intensive breaking waves in the surfzone, it was claimed in Kobayashi and Raichle (1994) that the number of events was satisfactorily reproduced in their model though admitting inherent uncertainties in the measured events such as noises by air bubble. The overtopping thickness had to be magnified by a factor of 5 to reasonably match with the measured one. However, the mean overtopping rate was at large underpredicted by the model. Kobayashi and Wurjanto (1989) and subsequently Dodd (1998) validated their NLSW models against the same data set of monochromatic wave overtopping of a sloping sea wall by Saville (1955). The seawall was situated on a shallow foreshore, which accommodated

Fig. 1. Model computational domain and hydraulic boundaries.

760

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

the conditions for wave breaking. Again, the average discharge was considerably underpredicted by both models, particularly for tests of high wave steepness (see details in Section 4.1). In this study, the numerical model based on the conventional NLSW equations as described earlier is applied to predict wave overtopping on a coastal barrier situated on a shallow foreshore (Tuan et al., 2006). In this case, wave breaking was intensive on the foreshore. The model results (see Section 4.2 for details) indicate that both the average discharge and the number of overtopping events are markedly underpredicted. In an interesting study on wave overtopping on a trapezoidal dike in the surfzone, Stansby and Feng (2004) reported that their NLSW model markedly underpredicted the mean water level in the surfzone compared to PIV measurements. The under-prediction appeared to be increasingly significant shorewards. The failure in reproducing the mean water level by the model was attributed to the inability of the NLSW model in describing the complex wave interactions in the surfzone, in which wave setup results. The computed water level in the surfzone in the absence of wave setup may partly explain the underestimation of wave overtopping by the aforementioned NLSW models. In addition, the friction coefficient can play a role. In general, it is rather straightforward to tune the friction coefficient in the numerical model in order to get the average discharge as measured because these are strongly correlated. However, it becomes complicated if we also consider the number of overtopping events, which is not sensitive to the change of the friction coefficient. To figure out which of these two factors, viz. the friction coefficient and the mean water level, will mostly affect the model results and thus help finding a solution to improve the shortcomings of the present conventional NLSW model, sensitivity analyses of the average overtopping discharge and the number of overtopping events by varying the bottom friction or the mean water level are carried out for a representative test selected from the aforementioned experiments by Tuan et al. (2006). The selected test was amongst the most intensive wave breaking tests (Iribarren number ξ = 0.55; generated

deep-water wave parameters: Hm0 = 0.163 m, Tp = 2.0 s; see Tuan et al., 2006), aiming to clearly demonstrate the effect of wave breaking on model predictions. The results are given in Fig. 2 for the average discharge and the number of events separately. For the effect of the bed friction, the friction coefficient fw is varied in the range between 0.0005 and 0.01 (for smooth slopes, see Section 4.2), while the mean water level is kept constant at 0.69 m as actually measured. For the effect of the water level, the mean water level is increased from 0.685 m to 0.71 m (2.5 cm), while the friction coefficient is fixed, viz. fw = 0.0035 based on the overall agreement of the average discharge between computations and measurements (see for details in Section 4.2). It should be noted that the increase of the water level here is of the same order of magnitude as the expected wave setup in the experiment. The results indicate that both the average overtopping discharge and the number of overtopping events are highly sensitive to the change of the mean water level (Fig. 2c and d), but much less to the change of the friction coefficient (Fig. 2a and b). It is emphasized that the sensitivity of the computed wave overtopping parameters to the change of the friction coefficient fw is not the same for all the tests. Therefore, the friction coefficient eventually chosen for the model must be based on overall agreement of the whole dataset (see details in Section 4.2). In the selected test case, wave breaking is so intensive so that the model seriously underestimates both the average overtopping discharge and the number of overtopping events even when fw approaches the value fw = 0.0005 for which the slope is considered as very smooth. In contrast, good agreement with the measurement would be obtained if the mean water level is increased by 1 cm, viz. from 0.69 m to 0.70 m. It should be noted that the measurement accuracy of the water level was in the range of ±1 mm during the experiment. Model computations, also for other tests of wave breaking, imply that the NLSW model underestimation is systematic and intrinsic. The deficit of the mean water level is in the range of the expected wave setup, which hints at that the model

Fig. 2. Effects of bottom friction and mean water level on average wave overtopping discharge q and number of overtopping events Novt. (a) and (b): Effects of bottom friction (c) and (d): effects of mean water level.

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

shortcoming might be due to the inability of modelling wave setup in the surfzone, i.e. wave breaking is not sufficiently modelled. This implication is in agreement with the observation by Stansby and Feng (2004) addressed earlier. During a storm surge, wave setup can be a very significant constituent of the total surge height (see e.g. Dean and Bender, 2006). Hence, the absence of wave setup in the mean water level can lead to a serious underestimate of wave overtopping. Overall, there are several factors discussed as follows that may be attributed to the aforementioned incapability of NLSW models. First, it is known that NLSW models are generally not able to describe wave frequency dispersion, thus producing the wave profile incorrectly, especially outside the surfzone (see Brocchini and Dodd, 2008). This consequently affects the model prediction of the characteristics of wave overtopping events accordingly. Second, an overtopping event is in fact the result of wave–structure interaction. However, this interaction could not be modelled properly in some cases such as steep (or near-vertical) slopes and slopes of complex geometry, where the underlying assumption about the hydrostatic pressure distribution of the NLSW flows is seriously violated. Third, in the case of breaking waves as quoted earlier and shown theoretically in Section 2.3, the effect of wave breaking is insufficiently modelled. This leads to an underestimate of the mean water level in the surfzone and consequently to a noticeable under-prediction of the mean overtopping rate and the number of overtopping events. It is generally realized that under fixed given hydraulic conditions the prediction of wave overtopping, owing to the non-dispersive behaviours of NLSW models, would considerably vary with the position of the imposed inflow boundary (Pullen and Allsop, 2003). The underestimation of both the average overtopping rate and the number of overtopping events becomes increasingly larger as the distance from the imposed inflow boundary to the dike toe increases. This is because the effect of frequency dispersion is larger seawards and becomes very significant outside the surfzone. Therefore, the first problem can effectively be resolved by imposing the inflow boundary conditions well within the inner surfzone and possibly closest to the dike toe; for instance, using a Boussinesq type model to transmit and generate wave conditions at this boundary. The second problem is inherent to a NLSW model and thus could not be resolved properly. However, the effect from a non-hydrostatic pressure on wave overtopping is considered to be small if the structure slope is not too steep and has a simple geometry like in the case of seadikes being considered. For slopes with complex structure configurations (e.g. with the presence of steep-faced walls), recent studies such as by Hu et al. (2000), Zhou et al. (2001) and Shiach et al. (2004) suggest that numerical and theoretical treatments can be applied so that NLSW models are still viable. Nevertheless, such findings are valid for the prediction of the average overtopping discharge only; no quantitative study is yet available for improvement of the prediction of individual overtopping events. In this paper, measures to deal with these two problems will not be further addressed as they are beyond the research scope. However, their influences on the model results will be discussed at some length later on. In Section 2.3, an approach to resolve the third problem is discussed in detail. 2.3. Source term related to energy dissipation by surface rollers In the surfzone, wave breaking is associated with the drastic motion of surface rollers, which dissipates the wave energy intensively. Therefore, the momentum flux balance through a breaking wave must be altered by the momentum flux carrying by the surface roller. However, this is not the case by the conventional NLSW equations. Preceding studies have indicated that wave breaking can more properly be modelled through incorporating an artificial breaking term related to the roller motion in the depth-averaged momentum equation. Svendsen and Madsen (1984) showed that the added breaking term attracts

761

energy from the roller motion and helps stabilize the wave front, which is an improvement over the conventional NLSW equations. Schäffer et al. (1993) pioneered in extending a Boussinesq model for breaking waves through including a roller-type convective term. Svendsen et al. (2000) confirmed experimentally that the breaking term is needed to maintain the momentum flux balance through a breaking wave and so wave breaking can more properly be modelled. In the literature, a number of approaches making use of the roller concept exist to describe wave breaking in both Boussinesq-type and NLSW models (e.g. Brocchini et al., 1992; Schäffer et al., 1993; Madsen et al., 1997). In the context of the later considered here, it is shown theoretically in the following that a roller term as such is indeed required to improve the model performance reduced by the underestimate of the mean water level as addressed in Section 2.2. This term is determined in connection with the energy dissipation by the surface roller. In the literature, it is widely accepted that the organized wave energy dissipated in breaking is first converted to the turbulent kinetic energy of surface rollers before being ultimately dissipated into small-scale turbulent motions. The surface roller is regarded as a mass of water, which rides on the underlying wave and stays stationary relatively to the wave crest (see e.g. Svendsen, 1984; Stive and De Vriend, 1987; Nairn et al., 1990; Tajima and Madsen, 2003, 2006). Hence, there must be a shear force between the surface roller and the underlying wave, which is apparently not zero in the surfzone as is commonly assumed by the free surface boundary condition of the NLSW equations. The (non-zero) surface shear force acts as an additional thrust toward the shore, gradually elevating the mean water level in the surfzone shorewards, viz. generation of (static) wave setup. The zero stress condition at the surface implies that wave setup is insufficiently described by the NLSW equations if surface rollers are present. Because wave overtopping is highly sensitive to the change of the mean water level, this may explain why NLSW models systematically under-predict both the average overtopping discharge and the number of overtopping events in the case of breaking waves as demonstrated in Section 2.2. It is delineated in the following that from the two dimensional RANS equations it is relatively straightforward to derive the horizontal momentum equation of NLSW flows, showing that an extra force balance term is needed if surface rollers are present. The two dimensional RANS equations, which for incompressible fluid are the continuity equation: ∂u ∂w + =0 ∂x ∂z

ð10Þ

and the momentum equation: ! ∂u ∂u ∂u 1 ∂p 1 ∂τxx 1 ∂τxz ∂2 u ∂2 u +u +w =− + + +ν 2 + 2 ρ ∂x ρ ∂x ρ ∂z ∂t ∂x ∂z ∂ x ∂ z ð11Þ in which u and w are the horizontal (x) and vertical (z) velocities, respectively, τxx and τzx are the turbulent (Reynolds) stresses, and p is the flow pressure. Neglecting vertical acceleration terms, the pressure distribution may be assumed hydrostatic so as to derive the NLSW equations: p = ρg ðh−zÞ

ð12Þ

Multiplying Eq. (10) by u and adding it to Eq. (11) yields: 2

2

2

∂u ∂u ∂uw 1 ∂p 1 ∂τxx 1 ∂τxz ∂ u ∂ u + =− + + + +ν 2 + 2 ρ ∂x ρ ∂x ρ ∂z ∂t ∂z ∂x ∂ x ∂ z

!

ð13Þ

762

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

After rearranging it follows:     ∂u ∂ p 1 ∂u ∂ 1 ∂u 2 = ð14Þ + u + − τxx −ν −uw + τzx + ν ρ ρ ρ ∂t ∂x ∂x ∂z ∂z where ν and νt are molecular and eddy viscosities, respectively. The Reynolds stresses can be determined according to the Boussinesq eddy-viscosity concept as follows: τxx τzx

R* = Rjp −

∂u = = ρνt ∂x ∂u = −ρu′ w′ = ρνt ∂z −ρu′ u′

ð15Þ

Integrating Eq. (11) over the flow cross-section A and neglecting the effect of molecular viscous stresses, the momentum equation is now transformed into: ∂U ∂ 1 A 2 + ð1 + βÞU + ∫ pdA− Txx ρA ρ ∂t ∂x

! =

1 ∂ ∫ τ dA ρ A ∂z zx

ð16Þ

in which the averaged quantities are defined according to:
2 • Energy correction factor: β =∫ u dA−1, β = 0 implies a uniform U A flow. 1 • Depth-averaged velocity: U = A ∫ udA A

with ν̄t • Depth-averaged normal stress: Txx = A1 ∫A τxx dA = ρνt ∂U ∂x being the depth-averaged eddy viscosity. Eq. (16) can be rewritten as: 2

∂Q ∂ Q + gI1P + ∂t ∂x A where

∂ ðβÞ≈0 ∂x

! =

From the above arguments, R⁎ resulting from the integration over the entire water depth including the roller may induce errors in the horizontal momentum balance. Therefore, it is prudent to treat the roller separately as a turbulence source regardless of its internal turbulent structure. To this end, R⁎ is integrated until the lower limit of the roller and then supplemented with the shear stress imposing by the roller:

  ∂ A 1 ∂ Txx + ∫ τzx dA ρ A ∂z ∂x ρ

ð17Þ

has been substituted and I1P is the hydrostatic presA

R* =

 1 ∂ 1 1
z = Z ∫ τ dA = τzx j z = 0s = B τzx j z = Zs −τzx j z = Zb ρ A ∂z zx ρ ρ

ð18Þ

where Zs and Zb stand for the elevations at the surface and the bottom, respectively. If the shear stress is assumed zero at the surface, R⁎ retains only the bed shear stress, i.e. the second term on the right hand side of Eq. (18). Further simplification by neglecting the contribution from the normal stress Txx, Eq. (17) exactly reduces to the ordinary horizontal momentum equation of NLSW flows as shown in Eqs. (1) and (2). It is worth noting that the integration in Eq. (18) is mathematically viable only if τzx is differentiable at any arbitrary elevation z over the flow depth domain. In other words, τzx must be first continuous and ∂τzx then must exist or be bounded at any z. ∂z However, this may not be the case over a vertical section in a breaking wave, where the roller slides against the underlying wave. In the mechanical sense, along the lower limit of the roller there exists a sliding stress against the flow underneath. This stress makes τzx discontinuous and thus non-differentiable along this boundary. Further, as the flow structure in the roller is violently turbulent with entrained air bubbles and vortices, it is uncertain whether τzx is continuous over the roller.

ð19Þ

in which the subscript (r) denotes the lower limit of the roller, τr is the shear stress imposed by the roller. The stress term Rjp arising from the roller reads: Rjp =

1 1A Bτ j = τ ρ zx z = Zs ρd r

ð20Þ

The roller is regarded as an independent source of turbulent energy, which transmits into the flow underneath. Mechanically, τr can thus be estimated through the work done by the roller. τr =

Dt D* = ir εj ur ur

ð21Þ

In which Dt (kg/s3 or W/m2) is the rate of production of turbulent energy density (per unit horizontal area) by the roller or the rate of loss of wave energy density by roller turbulence; Dir⁎ is the rate of dissipation of roller energy density, which generally lags behind Dt; ur is the relative (Lagrangian) velocity between the surface roller and the underlying wave; εj is an energy efficiency factor, for which εj = 1.0 implies a local equilibrium. As the wave crest and thus the roller both travel at the speed of the wave celerity in the Eulerian frame of reference of the computational domain, Eq. (21) can alternatively be expressed as follows (see also Stive and De Vriend, 1987):

sure force term I1P = ∫(h − z)dA. The second term on the right hand side of Eq. (17), hereafter designated R⁎, accounts for the effects of the turbulent shear stress. Without the presence of discontinuities (e.g. surface rollers) in the flow, this term can be quantified by integrating over the flow depth as follows:

1A 1A ðτ −τb Þ τ = ρd b ρd r

τr =

Dir c

ð22Þ

where Dir is the average rate of dissipation of roller energy density (hereafter referred to as the roller energy dissipation for brevity), c is the wave celerity. Substituting Eq. (22) into Eq. (19) and then amending to the conservative form as in Eq. (2), an additional term emerges in the source term vector as shown below: 0 1 0 ! →
→  → B
C S U; x = @ ∂ ∂U A + gA Sbx −Sf + Sr νt A ∂x ∂x

ð23Þ

where Sr is named as the roller slope term and is related to the roller energy dissipation Dir via: Sr = βr

Dir ρgdc

ð24Þ

where βr is the roller model coefficient of the order of 1. The first term in brackets on the right hand side of Eq. (23) arises from the horizontal velocity acceleration, which is negligibly small compared with other terms. The source term coefficient β̂ in Eq. (9) now includes Sr term and reads:
 ˆ = gΔx Aˆ Sˆbx − Sˆf + Sˆr β 2 cˆ

ð25Þ

where Ŝr is the Roe average of the roller slope term at a computational cell.

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

By definition, the roller slope term Ŝr is proportional to the amount of energy dissipation by rollers and thus becomes increasingly larger in the surfzone as wave propagates shorewards, creating a slope of the mean water surface like wave setup. In the case of non-breaking waves, this term automatically reduces to zero, and therefore has no influence on the model prediction. As an illustrative example of the influence of the roller term, time series of surface elevations computed with and without the roller term, for the same test in Section 2.2, are compared in Fig. 3 for three locations inside the surfzone (see also Section 2.4 for the determination of the roller term and Section 4.2 for details of the computation). It follows that with the roller term incorporated in the NLSW equations the water mean level in the surfzone is implicitly supplemented with wave setup and realistically elevates shorewards, which consequently improves model predictions of wave overtopping. Compared to the case without wave setup, variations of surface elevation are considerably changed, especially for locations well inside the surfzone. 2.4. Determination of roller dissipation Dir The roller energy dissipation Dir is the only unknown in Eq. (24) for determining the roller slope term Sr. There are two distinct approaches to implement this dissipation term into a numerical model. One is to apply Dir locally at the breaking wave front (see e.g. Schäffer et al., 1993 for an approach to locally include Dir as a convective term in a Boussinesq-type model for breaking waves). The value of Dir is instantaneously determined and the computation is therefore iterative since both bore height and location are not known in advance at a computational level.

763

Alternatively, it is assumed that the influence of the roller energy dissipation on wave overtopping follows a more static behaviour similar to wave setup. Therefore, Dir can be evaluated in a timeaveraged sense and can suitably be determined through a parametric (wave-averaged) wave energy decay model (termed as ENDEC hereafter), which solves energy balance equations expressing the energy exchange between the roller and the organized wave motion. For the present model we advocate the use of the latter approach because of its consistency with the problem of modelling the mean water level in the surfzone addressed in Section 2.2. To this end, the selection of a suitable ENDEC model is discussed below. The appropriateness of the selected strategy for incorporating the roller effect is further discussed in Section 3. The development of a class of parametric models initiated from the pioneering work of Battjes and Janssen (1978). In this work, the authors described the complex hydrodynamics of random wave breaking through two basic simple elements: a bore-type expression for the energy dissipation of breaking waves and an idealized distribution of the breaking wave heights. Since then the original model has been extended and validated against a vast number of datasets from both laboratory and field measurements, showing that the model is reliable for prediction of wave parameters including wave setup across the surfzone (see Battjes and Janssen, 2008 for a review on the history of developments of the model). Although several refinements and modifications have been made to the both model elements over the past decades, Battjes and Janssen (2008) pointed out that the revised versions produce only marginal improvements over the original one for most conditions. Perhaps the most important refinement so far is the incorporation of surface roller dynamics to improve the model

Fig. 3. Surface elevation variations with time computed with roller term (full line) and without roller term (dash-dot line): (a) at inflow boundary depth d = 0.16 m (b) at inner surfzone depth d = 0.05 m (c) in runup zone 0.015 m above still water level (horizontal dash line at 0.685 m).

764

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

prediction of the cross-shore profile of setup (or mean water level) (see e.g. Stive and De Vriend, 1994; Bosboom et al., 2000; Apotsos et al., 2007). Tajima and Madsen (2003, 2006) developed a similar parametric model, which is also governed by simple wave energy equations based on linear wave theory. However, the model deviates from the class of Battjes and Janssen's (1978) model in that the concept of equivalent linear wave characterising actual non-linear wave is used. Further, the energy dissipation due to wave breaking is determined according to the concept of stable wave energy instead of bore analogy (Dally et al., 1985). The evolution and dissipation of the roller energy across the surfzone are also implemented into the model through an energy balance equation in which the roller dissipation is determined in consistent with the dissipation due to breaking waves. Model predictions of the across-shore wave parameters compared well with measurements from a number of existing laboratory test cases. Again, the prediction of the across-shore setup profile was found to be considerably improved with the inclusion of the roller dynamics. Gonzales-Rodriguez (2006) carried out cross-comparisons between these two types of parametric wave model and found that their predictive capacity is very similar and hence either of these can be used for the present model. We adopt herein the model class of Battjes and Janssen (1978) because of its long history development as addressed in Battjes and Janssen (2008). The model elements are essentially based on the original model by Battjes and Janssen (1978) however with the (complete) Rayleigh distribution for breaking wave heights in the surfzone (see Janssen and Battjes, 2007). This modification, as stated in Janssen and Battjes (2007), is necessary to resolve the model unrealistic solution near the water line. In order to improve the model prediction of the cross-shore wave set-up profile, which is essential for the present research, the roller dynamics is incorporated. The roller energy and dissipation are determined according to Nairn (1990). For the sake of clarity, some basic formulations of this latest version of Battjes and Janssen's (1978) model are shown below: The wave energy balance equation reads:  ∂
ECg cosθ = −Dw −Df ∂x

ð26Þ

where Cg is the wave group velocity; θ is the wave incident angle (θ = 0 if wave propagates perpendicularly to the dike crest); Dw is the energy dissipation due to wave breaking; Df is the dissipation due to bottom friction (negligibly small compared to Dw); E is the wave energy defined according to linear wave theory. The energy dissipation Dw follows an explicit formulation by Janssen and Battjes (2007), assuming the complete Rayleigh distribution for wave heights in the surfzone: Dw =

pffiffiffi    3 
 3 π H 4 3 3 2 αfp ρg rms 1 + pffiffiffi R + R exp −R −erf ðRÞ 16 d 2 3 π ð27Þ

in which Hrms is the root-mean-square wave height; α is a model coefficient of the order of 1.0; fp = 1 / Tp is the peak (representative) frequency; erf represents the error function; R is the relative breaking depth defined as follows: R=

Hb γd = Hrms Hrms

ð28Þ

γ=

Hb = 0:39 + 0:56 tanh ð33S0 Þ d

ð29Þ

where d is the water depth; Hb is the breaker height; γ is the breaker index, which is a function of the deep water wave steepness S0.

The roller energy balance equation describing the evolution and dissipation of the roller energy in the surfzone reads (see Stive and De Vriend, 1994): ∂ ð2Er c cosθÞ = Dw −Dir ∂x

ð30Þ

where Er is the roller kinetic energy. It was assumed in Tajima and Madsen (2003, 2006) that only the loss of potential energy of broken waves, i.e. half of the total wave energy loss according to linear wave theory, contributes to the surface roller energy. Hence, in this special case the roller energy balance equation proposed by Tajima and Madsen (2003, 2006) is indeed identical to Eq. (30). The roller energy dissipation Dir is determined according to (Nairn, 1990): Dir = 2β0 g

Er c

ð31Þ

with β0 (=0.05–0.10) is the slope of the breaking wave face of random waves. The closed system of Eqs. (26)–(31) allows us to compute Dir at every computational cell. It is emphasized that the coupled ENDEC computation must begin well prior to the breaking point so the roller energy dissipation Dir can properly be determined. On the opposite, the inflow wave boundary for wave overtopping computations should be imposed after the breaking point and at the possibly closest location to the dike toe in order to minimize the effect from frequency dispersion not described the NLSW equations. 3. Discussion of roller term contribution The presence of the roller effect in the system of model equations has been established in previous sections. It remains to clarify the appropriateness of the inclusion of the roller term in the present context of modelling wave overtopping. In this section this problem is addressed regarding two main aspects: the magnitude and contribution of the roller dissipation term and its presence consistent with the present model. Regarding the model theory, the present model makes use of the flux-conservative form of the NLSW equations, which, in fact, satisfies jump conditions naturally (see Toro, 2001; Brocchini and Dodd, 2008). Therefore, the model is a shock-capturing one and intrinsically able to describe the process of wave breaking, viz. the momentum transferred from wave breaking to the water column is implicitly accounted for in the depth-averaged momentum equation through jump conditions. Despite the fact that wave breaking is naturally described through the numerical presentation of shock, physically, it is seen in Section 2.3 that when wave breaking with surface rollers an extra dissipation term is required in the momentum equation to account for the momentum exchange between surface rollers and the organized wave motion. An obvious issue to be discussed below is how we should introduce this term whether by its instantaneous values or time-averaged ones. In general, time- or phase- averaged quantities are often used in engineering applications instead of instantaneous ones. This is because, in a time- or phase-averaged sense, a sophisticated parameter can sometimes better be investigated and quantified whereas it is complicated and even not feasible at present to reliably measure and describe its instantaneous values. On this basis, Madsen et al. (1997) investigated the role of key parameters including roller terms of their Boussinesq-type (time-domain) model through averaging these quantities over a wave cycle. Regarding the roller effect, the analysis physically reveals that the contribution from the time-averaged roller dissipation term is mainly responsible for the cross-shore setup profile.

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

This is in line with observations from preceding studies, e.g. Stive and De Vriend (1994), Bosboom et al. (2000) and Apotsos et al. (2007), that inclusion of the roller dissipation considerably improve estimates of the wave setup profile across the surfzone. Further, an important finding from the work of Madsen et al. (1997) is that the time-averaged roller energy dissipation computed from their Boussinesq model is indeed in good agreement with that determined explicitly by the expression of Deigaard (1989) or Dally and Brown (1995). It should be noted that the expressions for the time-averaged roller dissipation according to Deigaard (1989), Dally and Brown (1995), and Nairn (1990) are essentially similar (see Eq. (31) and also Gonzales-Rodriguez, 2006), whereas the instantaneous roller dissipation in Madsen et al.'s (1997) model is according to Schäffer et al. (1993). From the above investigations two conclusions can be drawn out: first, the roller dissipation mainly contributes to create wave setup in the surfzone; second, the time-averaged roller dissipation determined from Boussinesq-type (using the roller concept) models and that by explicit formulations are of the same order of magnitude. Therefore, instead of applying the instantaneous roller dissipation term, it is reasonable and more efficient in the present context to consider the roller effect through wave setup so that the time-averaged roller dissipation term can be applied. In other words, we do not consider the roller effect on wave overtopping momentarily but on the scale of several wave cycles, i.e. through wave setup. The analysis in Section 2.2 of the model predictions sensitive to the change of the mean water level in comparison with existing laboratory data strongly supports this assumption, viz. lack of wave setup in the mean water level appears to be responsible for the poor model performance. It should be noted that the consideration of the roller dissipation term in a time-averaged sense in the present model is in essence similar to the incorporation of an independent equation describing the variation of the mean water level in the system of model equations. In this sense, the time-averaged roller contribution is consistent with the present time domain model. Last but not least, through incorporating the time-averaged roller dissipation term, besides the model efficiency, it is advantageous that calibration of the roller term in the present model comes down to calibration of the dissipation and evolution of surface rollers, which has already been resolved at length within the scope of the ENDECclass model itself (see Battjes and Janssen, 2008).

765

4. Model validation The developed NLSW model described in previous sections is validated herein against four laboratory datasets. The first two datasets of intensive breaking of regular and irregular waves (Sections 4.1 and 4.2) are used to see how well the supplemented roller term can help improving the model prediction of wave overtopping in the case of breaking waves. The last two datasets (non-breaking waves) are aimed at testing the model performance under more complex but realistic hydraulic conditions of combined surge overflow and wave overtopping (Section 4.3), and under a real grassed-dike slope (Section 4.4). 4.1. Monochromatic wave overtopping on a smooth impermeable seawall The well-known Saville (1955) small-scale data (Beach Erosion Board tests) of monochromatic wave overtopping of a smooth impermeable seawall are used here for validation of the present model supplemented with the roller term. The model results by Kobayashi and Wurjanto (1989) and by Dodd (1998), respectively referred to as KW and OTT models, validated against the same dataset are also presented for inter-comparison. The seawall in the experiment of Saville (1955) respectively had seaward slopes of 1:3 and 1:1.5 and was situated on a 1:10 foreshore. Monochromatic waves of various heights and periods were generated in a burst. The mean overtopping rate q for each test was determined over the stable portion of the incident wave train excluding reflected waves. Table 1 summarizes the parameters in prototype scale (upscaled by a factor of 17; see Dodd, 1998) of the 20 representative tests selected by Kobayashi and Wurjanto (1989). These parameters are kept the same here as were used by Kobayashi and Wurjanto (1989) and by Dodd (1998). In Table 1, Rc and hs are the dike crest freeboard and the water depth at the toe, respectively; Dt is the water depth determining the imposed location of the inflow boundary. The wave height Ht at the inflow boundary was determined from the deep water wave height H0 (at the wave paddle) using linear shoaling theory. Note that no measured wave profiles were available, instead the incident wave profile was specified using Cnoidal theory for Ur ≥ 26 (Ur is the Ursell number) and Stokes 2nd order theory for Ur b 26. This means that the actual wave profiles in the experiment might not correctly be reproduced in the model.

Table 1 Overtopping test parameters of Saville (1955) and model predictions. Test

H0 (m)

T (s)

Ur (m/s)

S0 (–)

Ht (m)

ξ (–)

Dt (m)

Rc (m)

hs (m)

Q measured (–)

Q KW (–)

Q OTT (–)

Qpres. (–)

Q pres. + bore (–)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

1.83 1.83 1.83 1.83 1.37 1.83 1.83 1.83 1.83 1.83 1.83 1.37 1.37 1.37 1.37 1.37 1.83 1.83 1.83 1.83

6.39 6.39 6.39 6.39 7.67 10.80 10.80 10.80 10.80 10.80 10.80 14.97 14.97 14.97 14.97 14.97 10.80 10.80 10.80 10.80

19.1 19.1 19.1 19.1 22.8 38.0 29.1 38.0 38.0 38.0 50.1 50.1 50.1 50.1 50.1 50.1 23.9 23.9 31.6 31.6

0.029 0.029 0.029 0.029 0.015 0.010 0.010 0.010 0.010 0.010 0.010 0.004 0.004 0.004 0.004 0.004 0.010 0.010 0.010 0.010

1.74 1.74 1.74 1.74 1.36 1.94 1.90 1.94 1.94 1.94 1.94 1.62 1.62 1.62 1.62 1.62 1.88 1.88 1.88 1.88

0.027 0.027 0.027 0.027 0.015 0.011 0.010 0.011 0.011 0.011 0.011 0.005 0.005 0.005 0.005 0.005 0.010 0.010 0.010 0.010

5.49 5.49 5.49 5.49 5.48 7.32 8.24 7.32 7.32 7.32 7.32 8.22 8.22 8.22 8.22 8.22 9.00 9.00 7.63 7.63

0.91 1.83 0.91 1.83 0.92 0.91 1.83 2.74 0.91 1.83 2.74 0.92 0.92 1.82 2.74 3.66 0.91 2.74 0.91 1.83

1.37 1.37 2.74 2.74 2.74 1.37 1.37 1.37 2.74 2.74 2.74 1.37 2.74 2.74 2.74 2.74 1.37 1.37 0.00 0.00

0.066 0.041 0.064 0.036 0.090 0.060 0.017 0.040 0.094 0.040 0.008 0.091 0.130 0.077 0.025 0.011 0.049 0.013 0.039 0.020

0.027 0.003 0.053 0.014 0.081 0.054 0.016 0.002 0.091 0.045 0.016 0.098 0.113 0.051 0.016 0.015 0.066 0.008 0.040 0.007

0.035 0.006 0.060 0.017 0.087 0.067 0.023 0.006 0.107 0.058 0.025 0.106 0.120 0.061 0.022 0.003 0.063 0.007 0.037 0.008

0.042 0.010 0.060 0.017 0.085 0.062 0.014 0.006 0.091 0.040 0.010 0.087 0.114 0.061 0.023 0.007 0.064 0.006 0.038 0.011

0.047 0.014 0.068 0.024 0.090 0.074 0.018 0.010 0.095 0.041 0.011 0.092 0.123 0.064 0.024 0.007 0.072 0.007 0.044 0.013

Seawall slope of 1:3 for tests 1 to 18, 1:1.5 for tests 19 and 20.

766

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

Fig. 4. Error of prediction of average discharge by KW, OTT, and present model.

For the sake of inter-comparison, the same friction coefficient fc = 0.05 as was used in KW and OTT is adopted here for computation. Fig. 4 plots the error of predictions of the dimensionless average qffiffiffiffiffiffiffiffiffi overtopping rate Q ( = q = gH03 ) as a function of the Iribarren number ξ (determined at the deep water boundary and with the beach slope 1:10) by KW, OTT, and the present model but without the roller term, respectively. The performance of the present model is comparable to that of OTT. All the three models appear to largely underpredict the average discharge (negative data points in Fig. 4). Moreover, the underestimation is inversely proportional to the value

of the Iribarren number ξ, or in other words, proportional to the degree of wave breaking in this case. The breaker depth, determined for each test using shallow shoaling theory with the breaker index γ estimated from Eq. (29), suggests that wave breaking does occur on the foreshore in all tests due to a depth-limited situation. Therefore, we may argue that this discharge underestimation is related to the insufficient description of wave breaking of NLSW models as shown theoretically in Section 2.3. Although the coupled ENDEC model described in Section 2.4 is for random waves, it is applied here to the case of monochromatic waves to see how well the model prediction of the average discharge can be improved. The average angle of the breaking wave face is taken as 14° (or in Eq. (31), β0 = tan14° = 0.25) as was used by Schäffer et al. (1993) for regular waves, while Eq. (29) is still assumed valid for the breaker index γ. In Fig. 5, model predictions of the dimensionless average discharge by the present model with and without the roller term in comparison with those by OTT and the measurements are shown. It follows that the present NLSW model with the roller term yields a reasonable improvement over those of conventional NLSW models. The modest improvement can be due to the fact that the roller term is determined on the basis of random breaking waves while the test case is monochromatic waves. Also, the large scattering of the model predictions (by both OTT and the present model) compared with the measurements is probably due to the difficulty in reproducing the actual wave time-series in the experiment as mentioned earlier. 4.2. Irregular wave overtopping on a smooth impermeable coastal barrier

Fig. 5. Average overtopping rates of monochromatic waves predicted with and without roller term.

Thirty five flume experiments of wave overtopping on a coastal barrier located on a shallow foreshore were carried out at the hydraulic laboratory of Delft University of Technology. The experiment summarized herein was described in detail in Tuan et al. (2006). The experimental setup is shown in Fig. 6. Two model barriers, made of smooth plywood, of different seaward (composite) slopes and crest heights were successively constructed on a gentle foreshore slope of 1/35. The tested waves were irregular (standard JONSWAP spectrum) with the zeroth-moment wave height Hm0 = 0.10–0.20 m and the

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

767

Fig. 6. Experimental setup of wave overtopping on a coastal barrier.

peak spectral period Tp = 1.5–2.7 s. The initial water depth in the flume was 70 cm in all tests. Under these experimental conditions intensive wave breaking occurred near the toe of the barrier and wave overtopping on the barrier was moderate to severe. The duration of each test was at least 1000Tp to assure the generation of the complete wave spectrum. The water level in the flume before and after each test was measured using a point gauge. Four capacitance wave gauges were used to measure the incident wave conditions, including those at the toe. A high sensitive wave gauge (WHM5) was accommodated in a well in the barrier crest to log the signal of overtopped waves, which was then used to count the number of overtopping events and the corresponding overtopping times. It was not aimed to measure the overtopping thickness in these small-scale experiments since this was considered unreliable (see Tuan et al., 2006). As addressed earlier, apart from the average overtopping discharge the characteristics of individual overtopping events are also important. It is referred to Tuan et al. (2006) for some newly proposed parameters to describe the characteristics of overtopping events. Here only two indicators are adopted: the number of overtopping events Novt and the relative overtopping time Fcd (≪1.0). The latter is defined as the ratio of the total overtopping time (flowing time of all overtopping events) to the storm (test) duration. A smaller Fcd implies that wave overtopping is more intermittent and vice versa. The measured wave time series at WHM4 (at the inner surfzone) together with the measured averaged water level in the flume for each test are used for the overtopping computation. The coupled ENDEC calculation begins at the deep section before the foreshore, i.e. at location of the 3-gauge array as shown in Fig. 6. The default value β0 = 0.10 for random waves is used. The friction coefficient fw generally has a significant influence on wave overtopping predictions of NLSW models. However, it is found in this particular case of intensive wave breaking that by tuning the friction coefficient alone does not resolve the problem of model underestimation. This is because the slope surface is already very smooth and the underestimation is too large (see Section 2.2). In general, the value of the friction coefficient must be determined experimentally according to overall agreement of wave overtopping data. In the literature, reference values are available for several types of “smooth” surface. Schüttrumpf and Oumeraci (2005) found experimentally that fw = 0.0058 for their smooth small-scale dike slopes made of wood fibreboard, which is much smaller than fw = 0.02 as suggested by Van Gent (1995) also for smooth slopes but constructed of concrete. The value fw = 0.01 for bare compacted clay slopes was recommended by Schüttrumpf et al. (2002) based on large-scale experimental results. The small-scale barrier slope considered herein was made of plywood so we may expect the friction coefficient fw is of the same order of magnitude as was suggested by Schüttrumpf and Oumeraci (2005) (fw = 0.0058). Based on the overall agreement of the mean overtopping rate between model predictions and measurements a somewhat smaller friction coefficient fw = 0.0035 is found, implying that plywood slopes are slightly smoother than those of fibreboard. Figs. 7 and 8(a and b) respectively depict predictions of the average discharge, the number of overtopping events and the relative over-

topping time, by both the model versions without and with the roller term, compared with the measurements. It follows that the conventional NLSW model underpredicts the average discharge, especially for the case of intensive breaking waves (cross points “×” with large deviations from the bisector in Fig. 7). Similarly, the number of overtopping events (Fig. 8a) and the relative overtopping time (Fig. 8b) are also underpredicted but more seriously, i.e. by 30% on average. Overall, the newly added roller term in the horizontal momentum equation appears to considerably improve the model predictions of the average overtopping discharge and the characteristics of overtopping events. Nevertheless, the number of overtopping events and the relative overtopping time are still slightly underpredicted by around 10%. A closer look at the model predictions of overtopping events as shown in Fig. 9, for the same test used in Section 2.2, indicates that some improvement on the spatial distribution of overtopping events is also achieved with the supplemented roller term. Interestingly, the computed and measured overtopping thicknesses are in the same range though it is doubtful about the quantitativeness of the measured values in the experiments.

4.3. Combined wave and surge overtopping on a levee Hughes (2008) carried out extensive small-scale experiments on combined (irregular) wave and (steady overflow) surge overtopping at earthen levees. The experimental setup is shown in Fig. 10. A model

Fig. 7. Measured and computed average wave overtopping discharge.

768

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

Fig. 8. Measured and computed wave overtopping events. (a) Number of wave overtopping events Novt (b) relative overtopping time Fcd.

levee was overtopped by steady surges superimposed with irregular waves under various (negative) freeboards. Seven pressure gauges (P1 to P7) were installed along the crest and the inner levee slope to measure the flow thickness. A LDV was also installed on the levee crest (at P2) to measure the time series of the horizontal flow velocity, which was then used to determine the

overtopping discharge at this location (as the product of velocity and overtopping thickness). The experimental programme was classified into 3 groups by the (negative) freeboard R = −1.2, −3.6, and −6.1 cm. In each group the tested wave parameters conforming to the JONSWAP spectrum were Hm0 = 3~11 cm and Tp = 1.2~2.9 s. In total, 27 experiments were carried

Fig. 9. Temporal distribution of wave overtopping events. (a) Measured events (b) computed events (conventional) (c) computed events (with roller term).

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

769

Fig. 10. Experimental setup of combined wave and surge overtopping on a levee (after Hughes, 2008).

out with various combinations of the freeboard and the wave parameters. The levee was initially overflowed by a steady surge and then followed by waves. Without loss of generality, only 9 of those experiments are considered herein (three experiments are picked up from each freeboard group). The measured surge level at P1 and the measured wave timeseries at the toe of the levee (WHM in Fig. 10) were used as input for computation. Since almost no wave breaking was expected on the foreshore of the levee the computation was performed with the conventional TVD schemes without coupling the roller term. The downstream boundary conditions were set to supercritical, like in the experiments. The friction coefficient fw = 0.0055, more or less the same as was suggested by Schüttrumpf and Oumeraci (2005), is found to give the best agreement with the measurements.

As an example, the time series of the flow depth and the horizontal velocity of test R131 (see Hughes, 2008) computed for the first 100 s at two locations, one on the levee crest (at P2) and one near the end of the inner slope (at P7), compared with the measurements are shown in Fig. 11. The computed parameters compare well with measurements, especially for the flow depth. Model predictions of the average discharge in the case of only steady current without waves (C) and in the case of steady current superimposed with waves (W + C) compared with the measurements from all 9 tests are shown in Fig. 12. It is noted that some uncertainties exist in the measurement of the horizontal velocity. These include signal dropouts (appeared as flat troughs in the velocity time series) during wave troughs as the flow depth drops below the laser beam of the LDV and that the measured velocity is not depth-averaged as computed. This may

Fig. 11. Combined wave and surge overtopping experiment R131: time series of measured versus computed flow parameters: (a) Flow depth at P2 (b) Horizontal velocity at P2 (c) Flow depth at P7.

770

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

4.4. Wave overtopping on a grassed seadike (GWK large-scale experiments)

Fig. 12. Combined wave and surge overtopping over a levee: measured versus computed average discharge.

explain why the agreement on the depth-averaged velocity and so on the average discharge is less satisfactory. The model results clearly demonstrate that the effect of waves on the average discharge decreases as the negative freeboard R (absolute value) increases. However, when the freeboard R is deeper than 3.6 cm this influence is already very small in the measured data whilst it is still considerable as predicted by the model. On the whole, the present model is very capable of simulating the flow on a dike under realistically complex hydraulic conditions of combined wave overtopping and surge overflow.

Under the framework of the European Eurograss and Floodsite projects, large scale (near prototype) experiments on dike breaching were carried out in the GWK wave flume at Hanover. The model dike was 5.8 m high and 5.0 m long (full flume width), and was composed of a sand core and a grassed-clay cover. The test programme was split into several stages, aiming to investigate wave forcing impact on the seaward slope, wave overtopping, and overtopping-induced breach initiation on the landward slope and breaching of the inhomogeneous dike. In this work, only overtopping experiments are considered, for which the experimental setup is shown in Fig. 13. Reference is made to Geisenhainer and Oumeraci (2008) for more details. Several wave gauge arrays were installed along the flume to determine the incident and reflected waves. The closest wave gauge (WHM in Fig.13b) was about 40 m from the dike toe. An overtopping chute together with a container mounted with a load cell was used to measure the overtopping discharge instantaneously (Fig. 13a). The width of the chute intake, adjustable to facilitate the overtopping measurement according to the expected overtopping amount in each test, was varied between 15 and 21 cm. The test programme consisted of 9 experiments with the generated wave parameters (TMA spectrum) at the wave board were the wave height Hm0 = 50~90 cm and the peak period Tp = 6~8 s. The flume water level was initially at 5.0 m in all the tests. The measured mean water level and the measured wave time series from WHM are used for the overtopping computation. The used friction coefficient fw is space-varying, viz. fw = 0.0055 for the (smooth) foreshore and 0.021 (or Chezy C = 30 m0.5/s) for the grass slope as suggested by Van der Meer (2007), respectively. Between these two sections, a transitional segment of about one wave length, over which fw increases linearly, is setup to avoid abrupt change of the bed friction. The effect of the roller term is negligibly small since there is almost no wave breaking on the dike foreshore in this case. In Fig. 14, the average discharges computed by the present model are shown in comparison with the measurements. The results for the number of overtopping events and the relative overtopping time are shown in Fig. 15, respectively.

Fig. 13. GWK experiments of wave overtopping on a grassed slope dike: (a) assembled chute and overtopping container (after Geisenhainer and Oumeraci, 2008) (b) experimental setup.

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

771

after propagating through the chute. The assembled chute and container clearly influence the capture of overtopping events. The overtopping time is significantly larger than the peak wave period for most measured events, which is unrealistic. It is very likely that the chute, whose slope was relatively flat, stretched out the overtopping events (i.e. prolonged the flowing time) and thus reduced the peak discharge accordingly. After all, it should be noted that in this case the inflow boundary (at WHM) is imposed well outside the surfzone. Therefore, one may expect some inaccuracies in the model predictions due to the inability of modelling frequency dispersion of the NLSW equations (see Section 2.2).

5. Conclusions

Fig. 14. GWK experiments: measured and computed average wave overtopping discharge.

The results show that both the average overtopping discharge and the number of overtopping events are fairly well predicted by the model. However, the model considerably underpredicts the relative overtopping time, i.e. the computed relative overtopping time is about one third of the measured one only (Fig. 15b). In other words, a computed event lasts considerably shorter than a measured one. Fig. 16 plots the computed and measured overtopping events together, showing marked differences in the order of magnitude of the peak discharge and the overtopping time of each event. It should be noted that the model computed overtopping events at a location on the dike crest whilst those from the experiments were measured at the container

A numerical model based on the nonlinear shallow water equations for wave overtopping on seadikes is presented. The model is extended with an extra source term related to the roller energy dissipation in the depth-averaged momentum equation to account for the additional effect resulting from the surface roller motion of breaking waves. The computation of wave overtopping thus requires the coupling with a parametric wave energy decay model to determine the dissipation and evolution of the surface roller energy in the surfzone. With the roller term, in essence, the computed mean water level is implicitly supplemented with wave setup, which consequently improves model predictions of wave overtopping. The results from model validation against existing laboratory datasets indicate that the supplemented roller term help improve model predictions of wave overtopping in the case of breaking waves considerably. The model also appears to be robust and capable of reliably simulating wave overtopping on seadikes under realistically complex hydraulic conditions. Because of the intrinsic drawbacks of the NLSW equations, it appears that NLSW models are generally more capable of predicting the mean overtopping rate than the characteristics of overtopping events. Besides, it is experienced that the capability of state-of-the-art techniques for capturing individual overtopping events is still far from sufficient, even in large-scale experiments, leaving behind uncertainties in measured

Fig. 15. GWK experiments: measured and computed wave overtopping events. (a) Number of overtopping events Novt (b) relative overtopping time Fcd.

772

T.Q. Tuan, H. Oumeraci / Coastal Engineering 57 (2010) 757–772

Fig. 16. GWK experiments: measured and computed temporal distribution of wave overtopping events.

data. This has greatly hampered our understanding of the wave overtopping hydrodynamics. Therefore, improvements of both numerical modelling capacity and laboratory techniques for more properly capturing individual overtopping events are urgently needed. Acknowledgements The research was funded by the Alexander von Humboldt (AvH) research foundation in collaboration with LWI-Braunschweig Technical University, Germany. The authors sincerely thank Professor Steven Hughes at Vicksburg Coastal and Hydraulics Laboratory (U.S. Army Engineer Research and Development Center) for providing valuable experimental data of combined wave and surge overtopping. Thanks are also due to Peter Geisenhainer and Piontkowitz for sharing the data from the GWK experiments. Constructive comments from anonymous reviewers helped improve the paper considerably. References Apotsos, A., Raubenheimer, B., Elgar, S., Guza, R.T., Smith, J.A., 2007. Effects of wave rollers and bottom stress on wave setup. J. Geophys. Res. 112. Battjes, J.A., Janssen, J.P.F.M., 1978. Energy loss and set-up due to breaking of random waves. Proc. 14th Int. Conf. Coastal Engineering. ASCE, pp. 466–480. Battjes, J.A., Janssen, T.T., 2008. Random wave breaking models: history and discussion. Proc. 31th Int. Conf. Coastal Engineering, Hamburg, Germany. Bosboom, J., Aarninkhof, S.G.J., Reniers, A.J.H.M., Roelvink, J.A. and Walstra, D.J.R., 2000. UNIBEST-TC 2.0 — Overview of model formulations. WL | Delft Hydraulics, Rep. H2305.42, Delft, the Netherlands. Brocchini, Maurizio, Dodd, Nicholas, 2008. Nonlinear shallow water equation modelling for coastal engineering. J. Waterw. Port Coast. Ocean Eng., ASCE 134 (2), 104–120. Brocchini, M., Drago, M., Iovenitti, L., 1992. The modelling of short waves in shallow waters: comparison of numerical models based on Boussinesq and Serre equations. Proc. 23th Coast. Engng. Conf., Venice, vol. 1. ASCE, Reston, VA, pp. 76–88. Burguete, J., Garcia-Navarro, P., 2001. Efficient construction of high-resolution TVD conservative schemes for equations with source terms: application to shallow water flows. Int. J. Numer. Meth. Fluids 37, 209–248. Dally, W.R., Brown, C.A., 1995. A modeling investigation of the breaking wave roller with application to cross-shore currents. J. Geophys. Res. 100 (C12), 24873–24883. Dally, W.R., Dean, R.G., Dalrymple, R.A., 1985. Wave height variation across beaches of arbitrary profile. J. Geophys. Res. 90 (C6), 11917–11927. Dean, R.G., Bender, C.J., 2006. Static wave setup with emphasis on damping effects by vegetation and bottom friction. Coast. Eng. 53, 149–156. Deigaard, R., 1989. Mathematical modelling of waves in the surf zone. Prog. Rep. 69. ISVA, Technical University, Lyngby, pp. 47–59. Dodd, N., 1998. A numerical model of wave run-up, overtopping and regeneration. J. Waterw. Port Coast. Ocean Eng., ASCE 124 (2), 73–81. EurOtop, 2007. Wave Overtopping of Sea Defences and Related Structures: Assessment Manual, Environment Agency UK/Expertise Netwerk Waterkeren NL/Kuratorium fur Forschung im Kusteningenieurswesen DE. Garcia-Navarro, P., Vazquez-Cendon, M.E., 2000. On numerical treatment of the source terms in the shallow water equations. Computer & Fluids, Elsevier 29, 951–979. Geisenhainer, P., Oumeraci, H., 2008. Large scale experiments on dike breaching and breach growth. Proc. 31th Int. Conf. Coastal Engineering. ASCE, Hamburg. Gonzales-Rodriguez, D., 2006. Modelling of nearshore hydrodynamics for sediment transport calculations. Master thesis, Massachusetts Institute of Technology, 133pp. Hu, K., Mingham, C.G., Causon, D.M., 2000. Numerical simulation of wave overtopping of coastal structures using the nonlinear shallow water equations. Coast. Eng. 41, 433–465. Hubbard, M.E., Dodd, N., 2002. A 2-D numerical model of wave run-up and overtopping. Coast. Eng. 47, 1–26.

Hughes, S.A., 2008. Combined wave and surge overtopping of levees: flow hydrodynamics and articulated concrete mat stability. Technical report ERDC/CHL TR-08-10, 184pp. Janssen, T.T., Battjes, J.A., 2007. A note on wave energy dissipation over steep beaches. Coast. Eng. 54, 711–716. Kobayashi, N., Raichle, A.W., 1994. Irregular wave overtopping of revetments in surf zones. J. Waterw. Port Coast. Ocean Eng., ASCE 120 (1), 56–73. Kobayashi, N., Wurjanto, A., 1989. Wave overtopping on coastal structures. J. Waterw. Port Coast. Ocean Eng., ASCE 115 (2), 235–251. Kobayashi, N., Wurjanto, A., 1992. Irregular wave setup and runup on beaches. J. Waterw. Port Coast. Ocean Eng., ASCE 118 (4), 368–386. Losada, I.J., Lara, J.L., Guanche, R., Gonzalez-Ondina, J.M., 2008. Numerical analysis of wave overtopping of rubble mound breakwaters. Coast. Eng. 55, 47–62. Madsen, P.A., SØrensen, O.R., Schäffer, H.A., 1997. Surf zone dynamics simulated by a Boussinesq type model. Part I. Model description and cross-shore motion of regular waves. Coast. Eng. 32, 255–287. Nairn, R.B., 1990. Prediction of cross-shore sediment transport and beach profile evolution. Doctoral thesis, Dept. Civil Eng., Imperial College, London. Nairn, R.B., Roelvink, J.A., Southgate, H.N., 1990. Transition zone width and implications for modelling surfzone hydrodynamics. Proc. 20th Int. Conf. Coastal Engineering. ASCE, pp. 68–82. Pullen, T., Allsop, W., 2003. Use of numerical models of wave overtopping: summary of current understanding. R&D Interim Guidance Note FD2410/GN1. Defra/Environment Agency. 7 pp. Reeve, D.E., Soliman, A., Lin, P.Z., 2008. A numerical study of combined overflow and wave overtopping over a smooth impermeable seawall. Coast. Eng. 55, 155–166. Roe, P.L., 1981. Approximate Riemann solvers: parameters and difference schemes. J. Comput. Phys. 43, 357–372. Saville, T., 1955. Laboratory data on wave run-up and overtopping on shore structures. TM-64, Beach Erosion Board, US Army Corps of Engineers, USA. Schäffer, H.A., Madsen, P.A., Deigaard, R., 1993. A Boussinesq model for waves breaking in shallow water. Coast. Eng. 20, 185–202. Schüttrumpf, H., Oumeraci, H., 2005. Layer thicknesses and velocities of wave overtopping flow at seadikes. Coast. Eng. 52, 473–495. Schüttrumpf, H., Möller, J., Oumeraci, H., 2002. Overtopping flow parameters on the inner slope of seadikes. Proc. 28th Int. Conf. Coastal Engineering, Cardiff. Shiach, J.B., Mingham, C.G., Ingram, D.M., Bruce, T., 2004. The applicability of the shallow water equations for modelling violent wave overtopping. Coast. Eng. 51, 1–15. Stansby, P.K., Feng, T., 2004. Surf zone wave overtopping a trapezoidal structure: 1-D modelling and PIV comparison. Coast. Eng. 51, 483–500. Stive, M.J.F., De Vriend, H.J., 1987. Quasi-3D nearshore current modelling: wave-induced secondary currents. Proc. Special Conf. Coastal Hydrodynamics. ASCE, New York. Stive, M.J.F., De Vriend, H.J., 1994. Shear stresses and mean flow in shoaling and breaking waves. Proc. Conf. Coastal Engineering, Kobe, Japan, pp. 594–605. Svendsen, I.A., 1984. Mass flux and undertow in a surf zone. Coast. Eng. 8, 303–329. Svendsen, I.A., Madsen, P.A., 1984. A turbulent bore on a beach. J. Fluid Mech. 148, 73–96. Svendsen, I.B.A., Veeramony, J., Bakunin, J., Kirby, J.T., 2000. The flow in weak turbulent hydraulic jumps. J. Fluid Mech. 418, 25–27. Tajima, Y., Madsen, O.S., 2003. Modelling near-shore waves and surface rollers. Proc. 2nd Int. Conf. on Asian and Pacific Coasts (on CD). World Scientific, Singapore. Tajima, Y., Madsen, O.S., 2006. Modelling near-shore waves, surface rollers, and undertow velocity profiles. J. Waterw. Port Coast. Ocean Eng., ASCE 132 (6), 429–438. TAW, 2002. Technical report wave run-up and wave overtopping at dikes. Technical Advisory Committee on Flood Defence, The Netherlands. 42 pp. Toro, E.F., 1997. Riemann Solvers and Numerical Methods for Fluid Dynamics: A Practical Introduction. Springer-Verlag, Berlin. 624 pp. Toro, E.F., 2001. Shock-Capturing Methods for Free-Surface Shallow Flows. Wiley, New York. 309 pp. Tuan, T.Q., Verhagen, H.J., Visser, P.J., Stive, M.J.F., 2006. Wave overwash at low-crested beach barriers. Coast. Eng. J., World Scientific and JSCE 48 (4), 371–393. Van der Meer, J.W., 2007. Workpackage 3: Development of Alternative OvertoppingResistant Sea Defences, Phase 3: Design, Construction, Calibration and Use of the Wave Overtopping Simulator. Rijkswaterstaat, Delft. Van Gent, M.R.A., 1995. Wave interaction with permeable coastal structures. Doctoral dissertation, Delft Hydraulics Press, 191 pp. Zhou, J.G., Causon, D.M., Mingham, C.G., Ingram, D.M., 2001. The surface gradient method for the treatment of source terms in the shallow water equations. J. Comput. Phys. 168, 1–25.