MODELLING OF WAVE OVERTOPPING FLOW OVER COMPLEX DIKE GEOMETRIES: CASE STUDY OF THE AFSLUITDIJK

Grass cover erosion by overtopping waves is one of the main failure mechanisms of dikes. Transitions in cover type and geometry can increase the hydraulic load and are therefore identified as vulnerable locations for grass cover erosion. Two models are applied to the inner slope of the Afsluitdijk in the Netherlands to show how transitions can be included in overtopping models. Firstly, the analytical grass-erosion model is used to simulate the erosion depth along the profile for a six-hour storm. The model results show that the erosion depth is maximal at the end of the two slopes in the profile. Secondly, the effect of transitions on the hydraulic load is computed with a detailed hydrodynamic model. The model results show that geometric transitions significantly influence the shear stress, the normal stress and the pressure. Four vulnerable locations for grass cover erosion are identified based on the model results that are related to slope changes along the profile. Furthermore, the model results show that the overtopping flow is mainly affected by geometric transitions, while no effect of roughness transitions on the modelled forces was observed.


INTRODUCTION
Wave overtopping is one of the main failure mechanisms of grass-covered dikes. The process works as follows: water waves run up the waterside slope, overtop the crest and subsequently accelerate on the landward slope. Although good grass covers on top of a clay layer have shown that they are able to withstand large overtopping volumes, locally the grass cover may erode due to large flow velocities and high amounts of turbulence ( Figure 1). A complex dike geometry is a dike that contains several transitions such as transitions in cover roughness, in geometry or an existing erosion hole. These transitions affect the flow of the overtopping wave and thereby the grass cover erosion. For example, roughness differences can lead to additional turbulence and jets can form at geometric transitions ( Figure 2). These transitions are vulnerable spots for grass cover failure as experiments have shown that grass erosion often starts at transitions (Steendam et al., 2014). Therefore, it is important to gain knowledge on the hydraulic forces resulting in erosion at these transitions. The additional erosion at transitions can be caused by an increase in hydraulic load on the dike cover, a decrease in the cover strength or a combination of both.
As a first step, detailed information on the hydraulic load at transitions is required. Two types of loading are described in models for grass-cover erosion: normal forces (Ponsioen et al., 2019) and shear forces (Dean et al., 2010;Van der Meer et al., 2010;Hoffmans, 2012). The normal force mechanisms describes the vertical pull of the overtopping flow on the grass cover leading to under pressures in the cover layer and lifting of the dike cover ( Figure 3). The overtopping flow also pulls horizontally on the grass cover at locations with high (shear) velocities and high amount turbulence. Therefore, the normal stress, the pressure, the shear stress and the flow velocity are necessary variables to describe the hydraulic load at transitions. U τ n τ s Figure 3: Schematization of the a wave overtopping the dike crest indicating the maximum flow velocity U, the normal stress τ n working normal to the dike surface and the shear stress τ s working parallel to the dike surface.
Current calculation methods for overtopping failure are often based on the maximum flow velocity in 1D erosion models across the dike, such as the cumulative overload method (Van der Meer et al., 2010) and the analytical grass-erosion model (Van Bergeijk et al., 2019a). Several formulas for the maximum flow velocity along the dike profile exist (Schüttrumpf and Oumeraci, 2005;Hughes et al., 2012;Van Damme, 2016;Van Bergeijk et al., 2019b). However, transitions not solely affect the flow velocity of the overtopping wave, but also the amount of turbulence in the overtopping wave. Therefore, these erosion models require an additional factor for the load at transitions (Van Hoven et al., 2013;Van Bergeijk et al., 2019a). The values of these factors are not known for all transitions and theoretical formulations for these factors can underestimate the additional load at transitions (Warmink et al., 2020a).
Computational fluid dynamics models can compute the flow velocity, pressure, shear stress and normal stress of the overtopping waves (Bomers et al., 2018;Van Bergeijk et al., 2020). Additionally, the turbulence is solved within the models and therefore the effect of transitions on the hydraulic load can be studied directly. Previous model studies have shown that geometric changes are vulnerable locations for grass cover erosion due to an increase in near-bed turbulence (Bomers et al., 2018) and high shear forces (Van Bergeijk et al., 2020). However, solely one type of transition has been studied with these models and the effect of multiple transitions along a dike profile is still unknown.
This study shows the application of two newly developed overtopping models to the inner slope of the Afsluitdijk in the Netherlands, which contains several transitions in geometry and in cover type. The goal of this study is to identify vulnerable location for grass cover erosion by overtopping waves based on the hydraulic load near transitions. First, the analytical Grass-Erosion Model (GEM) is applied to the Afsluitdijk to compute the erosion depth along the dike profile. Two formulations for the turbulence parameter are studied to show the effect of an additional load factor at transitions in a flow-velocity based erosion model. Next, the recently developed 2DV model for the wave overtopping flow over a grass-covered dike (Van Bergeijk et al., 2020) is used to compute the flow velocity, pressure, shear and normal stress along the dike profile for four overtopping volumes. The hydraulic load along the dike profile is analyzed to gain insights in the effect of transitions on the overtopping flow and the hydraulic load near transitions. Four vulnerable locations for cover erosion are identified in this case study. In the end, the model results are compared to other studies of the inner slope of the Afsluitdijk.

CASE STUDY: AFSLUITDIJK
The Afsluitdijk is a dam between the Lake IJssel and the Wadden Sea in the Netherlands (Figure 4). The Afsluitdijk is currently reinforced to meet the new Dutch safety standards. The outer slope on the Wadden Sea side is protected against breaking waves by a combination of blocks, hard revetments and a berm (Klein Breteler et al., 2019). In the central section of the Afsluitdijk, the crest and inner slope on the Lake IJssel side are covered in grass, but also contain two berms with a biking path and two roads for transportation ( Figure 4b). This results in multiple transitions on the inner slope, both in cover type as well as in geometry.
No measurements of the exact dike profile are available, therefore the dimensions of the dike profile are based on Hoffmans (2014) with a smooth grass-covered crest and slope (Figure 5a and 6e). The biking path and roads are modelled by changing the roughness factor for the dike cover without height difference between the grass cover and the asphalt.
Experiments with the Wave Overtopping Simulator were carried out at dike section at the eastern end of the Afsluitdijk to study the strength of the grass cover (Bakker et al., 2009). Two dike sections were tested: (1) regular grass-covered dike profile, and (2) a grass-covered crest and slope with a road 3 m from the inner toe. At the first test section, erosion started at the inner toe for an overtopping discharge of 10 l/s/m ( Figure 1) and on the slope at 50 l/s/m. The inner toe at the second section was again damaged at 10 l/s/m and the road at 30 l/s/m.
Translating these results to the central section of the Afsluitdijk, we can expect that the transition from slope to berm and the transitions from grass to asphalt at the biking path and roads are vulnerable locations for erosion.

MODEL SET-UP
Two types of models for the overtopping flow over the crest and inner slope are used in this study. The analytical grass-erosion model is a fast model that calculates the flow velocity and erosion depth of every overtopping wave during a storm. This model provides insights in the amount of erosion during a storm and the location where most erosion occurs. Next, a detailed hydrodynamic model in OpenFOAM is used to simulate the hydraulic forces working on the dike cover and identify locations where the hydraulic load is high. Both models can be used to identify vulnerable locations for grass cover erosion along a complex dike geometry.

Analytical grass-erosion model
The analytical Grass-Erosion Model (GEM) calculates the erosion depth along the dike profile for all waves overtopping during a storm Van Bergeijk et al. (2019a). In this model, the load is described by the maximum flow velocity U and the turbulence parameter ω. The strength of the grass cover is described by the threshold flow velocity U t .
The erosion depth along the dike profile is computed for a storm with a duration of 6 h and an average overtopping discharge of 5 l/s/m. The threshold flow velocity is set to 19.2 m/s corresponding to a critical flow velocity of 8 m/s that was determined form the overtopping tests (Van Hoven, 2014).
The flow velocity along the inner slope of the Afsluitdijk is computed with the analytical formulas of Van Bergeijk et al. (2019b) and depends on the cover roughness and the slope angle. Transitions in cover type are simulated by adapting the friction factor f for the bed roughness from 0.01 for grass to 0.02 for asphalt at the biking path and roads. Moreover, transitions can influence the amount of turbulence, and thereby increase the load on the dike cover ( Figure 2). The erosion depth is computed for two turbulence formulation. First, a constant turbulence parameter ω of 2.0 is used. In this case, the effect of transitions on the hydraulic load is solely modelled using in the formulas of the flow velocity. Next, the effect of transitions are also included in the turbulence parameter using the formulas of Hoffmans (2012). In this case, the turbulence parameter depends on the slope angle and the friction factor, and therefore changes along the profile. The erosion depth along the dike profile is computed for both turbulence formulations and the locations where the erosion depth is maximal are identified.

Detailed hydrodynamic model
The 2DV model is developed in the open-source software OpenFOAM and simulates the flow from the start of the dike crest to the end of the inner slope (Van Bergeijk et al., 2020). The Nikuradse roughness height K s of grass was calibrated to 8 mm and the roughness of asphalt is set to 10.1 mm, corresponding to a friction factor f of 0.01 and 0.02 for a layer thickness of 20 cm, respectively. The horizontal grid size is 5 cm and the vertical grid size is 3 cm. The required boundary conditions at the start of the crest are the depth-averaged flow velocity and layer thickness of the overtopping wave as function of time which are generated from the overtopping volume (Van Bergeijk et al., 2020). The k − ω SST turbulence model is used to simulate the turbulence in the flow.
The hydraulic load along the dike profile is computed for four overtopping volumes V: 1, 2, 3 and 4 m 3 /m. The model output includes the flow velocity in the top layer, the pressure on the dike cover, the shear stress on the dike cover and the normal stress on the dike cover, which are all maximum at the wave front (Van Bergeijk et al., 2020) (Figure 3). The maximum flow velocity, maximum normal stress, maximum shear stress and maximum pressure are computed at every location along the dike profile for four individual volumes and the locations where the hydraulic variables are maximal are identified.

Flow velocity and erosion depth along the dike profile
The flow velocity increases along the slope and decreases along horizontal parts of the dike profile in the analytical grass-erosion model (Figure 5b). This results in a high load at the end of the upper slope (x = 11 m) and the end of the middle slope (x = 17 m). Therefore, the erosion depth is maximal at these locations ( Figure 5c). The formulation of Hoffmans (2012) results in a higher turbulence parameter ω on the slopes compared to the constant value of 2. On the horizontal parts, the turbulence parameter depends on the flow velocity and is higher than 2 at the start of the berms (x = 11 and 18 m), but lower than the constant value at the end of the berms at x = 15.5 and 27-31 m. This can be seen in the dashed line of the erosion depth that is larger when the the turbulence parameter is higher than the constant value, while the erosion depth is smaller when the turbulence parameter is smaller than the constant value.
The turbulence parameter does not influence the location of most erosion, which seems to be determined by the flow velocity in the GEM model. However, the turbulence model influences the amount of erosion. Since failure in the Dutch safety standards is defined as an exceedance of 20 cm erosion depth ('t Hart et al., 2016), the turbulence parameter formulation can determine if the cover does or does not fail for wave overtopping. Therefore, it is important to improve the formulations of the turbulence parameter and other additional factors in flow velocity based erosion models.

Forces on the dike cover
The results of the detailed hydrodynamic model show that the four hydraulic variables increase with the overtopping volume V (Figure 6). Similar to the GEM results, the maximum flow velocity U increases over the slopes (x = 3 − 11 m, 16-17 m and 41-45 m) and decreases along horizontal parts of the dike profile ( Figure 6a). The flow velocity is maximal at the end of the upper slope, thus most erosion is expected at this location using a flow velocity based erosion model.
The shear stress τ s is small on the dike crest and increases along the upper slope (Figure 6b). At the transition from the upper slope to the first berm (x = 11 m), the shear stress shows a small decrease followed by a sharp increase. The same pattern is observed on the middle slope and the second berm where a sharp decrease in followed by a sharp increase in shear stress (x = 16 − 18 m). The steepness of these peaks is larger for the two larger overtopping volumes (V = 3 − 4 m 3 /m) than for the smaller overtopping volumes (V = 1 − 2 m 3 /m). The shear stress is maximal around x = 18 m at the beginning of the second berm followed by a decrease along the profile where the shear stress becomes negligible small at the end of the profile for the lower two volumes. The results show that the modelled shear stress is not affected by the transitions in cover type between grass and asphalt.
The normal stress τ n is negligible small on horizontal parts of the profile and high on the upper and middle slope (Figure 6c). The location of the maximum normal stress depends on the overtopping volume: the upper slope for 1 m 3 /m and 3 m 3 /m, the middle slope for 2 m 3 /m and the start of the second berm for 4 m 3 /m. These locations are all related to a slope (change) which shows that the geometry has a large effect on the normal stress.
The pressure increases significantly at the slope changes around 3, 11, 16 and 17 m with a factor 3-6 depending on the overtopping volume (Figure 6d). Only for the smallest overtopping volume of 1 m 3 /m, the pressure decreases away from the middle slope (x = 18 m).

VULNERABLE LOCATIONS FOR GRASS COVER FAILURE
The model results show four locations along the profile that are vulnerable for grass cover failure due to a high hydraulic load of the overtopping wave (Figure 7): 1. Crest to upper slope: significant increase in the maximum pressure.
2. Upper slope: high flow velocity and shear stress. The normal stress is maximal for V = 1 m 3 /m and 3 m 3 /m.
3. Upper slope to berm: significant increase in pressure and the location where pressure and flow velocity are maximal as well as the normal stress for V = 2 m 3 /m. The GEM model predicts a large erosion depth at this location.
4. Middle berm: significant increase in pressure and the location where the shear stress is maximal. The normal stress is maximal for V = 4 m 3 /m. This location is especially vulnerable because of the two geometric transitions close to each other at a cross-dike distance of 16 and 17 m. Additionally, the GEM model predicts a large erosion depth at this location.
It is interesting to notice that these locations are all related to geometry such as a slope change or a slope itself. Moreover, we noticed that the hydraulic variables in the OpenFOAM are not significantly influenced by the transitions in cover type between grass and asphalt ( Figure 6). This could be related to the insensitivity of the numerical model to the roughness height (Van Bergeijk et al., 2020), but another more likely reason might be that the roughness difference between grass and asphalt is quite small as can be seen in the small change of the flow velocity in the GEM model when the friction factor changes at the grass -asphalt transitions. Therefore, it might be that roughness transitions only have a small contribution to turbulence generation and thereby a small contribution to the hydraulic load. (1) crest to upper slope, (2) upper slope, (3) upper slope to first berm, (4) middle slope between the two berms.

COMPARISON TO OTHER STUDIES
Wave overtopping tests on the Afsluitdijk have been performed in 2009 at a section without the berm with the cycling path (Bakker et al., 2009). These tests showed that the inner toe -comparable to location 3 in Figure 7 -is the weakest point along the profile. Moreover, significant erosion occurred at the start of the road (x = 19 m in Figure 6e) and the slope was the strongest location along the profile. Our model results show an increase in the hydraulic load at the inner toe which can explain the erosion at the inner toe. This is the location where the hydraulic variables are at their maximum, thus, most erosion is expected. However, the hydraulic load does not increase at the transition on grass to the road at the cross-dike distance x = 19 m. The measured erosion can be caused by a reduction of the cover strength at this location, or the difference in the profiles results in smaller forces in our modelled case due to the additional berm with biking path compared to the test section without this berm.
Furthermore, we assumed no height difference between the grass cover and asphalt roads which could be another reason why the hydraulic load does not change at the transitions in cover type. As seen in this study, geometrical transitions have a large influence on the hydraulic load, and therefore a height difference near a road could lead to additional load. This could be an explanation why little erosion was observed at overtopping experiments with a road on top of a grass-covered dike in Belgium (Peeters et al., 2012) with a smooth transition between the grass and asphalt, contrary to the same type of experiment on a dike in the Netherlands with a large height difference at this transition leading to a large amount of erosion (Van Hoven et al., 2013). Further research into the effect of height differences on the hydraulic load is necessary to confirm this hypothesis.
In this study, a six-hour storm was simulated in the analytical grass-erosion model using constant hydraulic conditions. This approach is commonly used to simulate wave overtopping on river dikes (Van Hoven, 2014). However, the storm conditions for a sea dike such as the Afsluitdijk depend not solely on the storm duration but also on the tides. Therefore, a storm event for a sea dike has a longer duration and the hydraulic conditions increase over time towards the peak conditions followed by decrease. Warmink et al. (2020b) incorporated the storm developed in the analytical grass-erosion model and applied the model to the section of the Afsluitdijk where the tests were performed. Warmink et al. (2020b) showed that a three-hour storm approach underestimates the hydraulic load compared to the storm development approach. Therefore, it is recommended for safety assessment purposes to include the storm development in the analytical grass-erosion model for sea dikes.

CONCLUSIONS
Complex dike geometries include several transitions in cover type and geometry that can influence the hydraulic load of the overtopping wave and the resulting grass cover erosion. Two overtopping models are applied to the inner slope of the Afsluitdijk in the Netherlands that contains several transitions. The model results are used to identify vulnerable locations for cover failure based on the hydraulic load at transitions.
First, the analytical grass-erosion model predicts most erosion at the end of the upper and middle slope. These locations are not influenced by the formulation of the turbulence parameter, but the erosion depth increases when the effect of the cover type and slope angle on the turbulence is included in the model.
Next, the detailed hydrodynamic model developed by Van Bergeijk et al. (2020) was used to gain insights in the forces working on the dike cover and to study the hydraulic load at transitions. The models results show showed that geometric transitions have a large influence on the hydraulic forces. The pressure and normal stress increase significantly at the slope changes where the pressure increases up to a factor 6. The shear stress increases on the slope and showed a large increase on the middle slope with two geometric transitions close to each other.
Both model results indicated four vulnerable locations for grass cover failure: (1) the transition from the crest to the upper slope, (2) the upper slope, (3) the transitions from the upper slope to the berm, and (4) the middle berm. These model results agree with the results of overtopping tests on the Afsluitdijk in 2009 that showed that the erosion started at the transition from the upper slope to the berm.
Changes in cover type have a minor influence on the forces in the numerical model and on the erosion depth in the analytical model. Larger roughness differences need to be studied to confirm that transitions in cover type do not affect the hydraulic forces.