MODELLING THE EFFECT OF COASTLINE ORIENTATION ON STORM ERTOSION AT THE SYLT ISLAND, NORTH SEA

The effect of coastline orientation on storm erosion was investigated using a hybrid numerical modelling approach (Delft3D and XBeach) on the convex-shape North Sea coast of Sylt. Storm impacts were simulated using different severity events occurring in a cluster and alone. Both models predicted reasonable wave dynamics (RMSE of wave spectral density ≤ 0.45 m/Hz). XBeach resulted in higher erosion above MSL than Delft3D. High severity events approached the coast from the SWNW sector (i.e., dominant storm approach on this coast). The convex-shape of the coastline provided however different storm landfall at individual locations along the coast. The coastal stretch from the apex to the North showed high vulnerability to storm erosion with the highest erosion at the apex, when the storm events approach from the northwesterly directions. This coastal stretch draws attention in the application of management strategies to withstand the storm erosion.


INTRODUCTION Background
Storm erosion along beach and dune systems mainly depends on 1) Severity of storm events, 2) Nearshore morphology and 3) Exposure of coastline towards the orientation of storm landfall. Storm severity indicates the erosive power of a storm event and can be estimated using a proxy storm power index, which is a function of storm wave heights and duration (Dissanayake et al., 2015b). The higher the index the larger the erosion during storm impacts. The presence of nearshore morphological features plays a dominant role against storm erosion, e.g., the ridge and runnel pattern of the macro-tidal environment, Sefton coast (Liverpool bay, UK), provides resistance to storm events approaching to the coast (Dissanayake et al., 2015a). The chronology of events in a storm cluster thus determines the amount of erosion. Smallegan et al. (2016) showed that the presence of a buried seawall increased wave attenuation and thus the resilience of the beach and dune system during the impacts of Hurricane Sandy at Bay Head, USA. Coastline orientation regulates the storm landfall at individual locations along coasts, which results in different erosion and sedimentation patterns. Guisado-Pintado and Jackson (2019) showed the variability of storm impacts along Five Finger Strand (NW, Irish coast) during two storm events with different severities. In addition to the synchronising of forcing, the coastline orientation influenced the erosion along this coast. Storm impacts on the Sefton coast showed that increased erosion at the apex and then sedimentation at the North and the South therein (Dissanayake et al., 2015c). In the present study, we focus on the impacts of coastline orientation on storm erosion along the convex-shape coastline of the Sylt island at the North Sea coast (Fig. 1). This coastline provides locally different orientations towards an approaching storm wave, and thus allows a systematic study of the effect of coastline orientation on storm erosion (Dissanayake et al., 2019).
A priori knowledge of potential storm erosion of a coastal system is of utmost importance to implement suitable management strategies against storm impacts (Vousdoukas et al., 2020). Along convex-shape coastlines, such knowledge is even more crucial because it is expected high alongshore variability of erosion compared with that of straight coastlines. Routine beach nourishments are carried out on the North Sea coast of Sylt, which intends to counteract the seasonal storm erosion (Blossier et al., 2017). Results of the present study might support to improve the nourishment quantities along the coast by identifying erosion volumes during different severity storm events.
Numerical modelling is increasingly used to investigate the storm impacted beach and dune erosion (e.g., Van Ormondt et al., 2020). In this study, we use two open-source morphodynamic models (Delft3D and XBeach) and compare storm erosion between models during three storm events, which have different level of severities (i.e., storm power index). Both models are simulated applying the default parameter setting for hydrodynamics. Bed sediment composition is implemented using a single sediment fraction (D50 = 0.2 mm).

Objective
The main objective is to investigate the spatial variability of storm impacted beach and dune erosion of a convex-shape coastline using two different models. Furthermore, model skills in hindcasting wave dynamics and morphodynamics are investigated based on the occurrence of different severity storm events.

Study area
The Sylt island is located in the North Frisian Wadden Sea (Fig. 1). The western coast of Sylt (North Sea coast) is about 35 km in length and its orientation from south to north varies about 20 o being the southern coast almost aligned with north. This study focused on the central coast, which encloses the stretch of the maximum change in coastline orientation. The nearshore morphology generally exhibits a double-barred profile with steep slope on the beach and regular rip channels (Blossier et al., 2017). Sediment is characterised by medium to coarse sand (Ahrendt, 2001). The dune system on Sylt reaches up to about 30 m in height. Beach nourishments are routinely carried out to combat the storm impacted beach and dune erosion (Blossier et al., 2017).
The semi-diurnal meso-tide has a maximum tidal range of about 2.5 m (at H in Fig. 1). The significant wave height (Hs) during storm events sporadically exceeds 6 m, while the mean Hs fluctuates around 1 m (at 13 m water depth: B in Fig. 1). The dominant wave direction is from NW, and the waves are mainly generated from the westerly wind from SW to NW (at LW in Fig. 1). Data from these observations were used to simulate the models. Additionally, the wave data at L (8 m water depth) were adopted to assess the model skills in predicting wave dynamics, which is the main driver for sediment transport and morphodynamics.

Models
A hybrid-approach of numerical modelling used two process-based morphodynamics models, Delft3D and XBeach, to simulate the storm impacts under different severity events, and to investigate the effect of coastline orientation.

Delft3D (D3D)
The D3D model suite was used in the depth-average approach (2DH) applying Hydrodynamics (FLOW), Sediment transport (SED) and Morphodynamics (MOR) modules together online-wave coupling with SWAN. SWAN is a 3 rd generation wave model, which estimates generation, propagation and dissipation of short-wave spectrum (Booij et al., 1999). In the flow module, the unsteady shallow water equations are solved by the alternating direction implicit method to compute the hydrodynamics (Stelling and Leendertse, 1991). Total sediment transport under combined waves and currents is estimated using the Soulsby -Van Rijn formulations (Soulsby, 1997). Morphodynamics are computed based on the conservation of sediment fluxes, which is multiplied by a morphological acceleration factor, morfac (Lesser et al., 2004).

XBeach (XB)
XB is a 2DH morphodynamic model, which estimates the beach and dune response to time-varying storm conditions. It includes the main physical processes of dune erosion by swash, collision, overwash, and inundation. XB solves coupled depth-average equations for wave propagation, flow, sediment transport and morphodynamics. The wave solver is based on the 2 nd generation HISWA model and computes the propagation of short wave averaged envelop, long-wave motion (Roelvink et al., 2009). Similar to D3D, hydrodynamics are computed using the shallow water equations. Avalanching is used to estimate dune erosion, and cross-shore transport depends on the onshore-component by wave skewness and asymmetry, and the offshore component by return flow. Similar to D3D, sediment transport is estimated using the Soulsby -Van Rijn formulations, and morphodynamics can be accelerated using a morfac.  Table 1. Both models were simulated using the default parameter settings for the numerical solver, and the real time storm impacted erosion was estimated without morphodynamic acceleration, morfac = 1.

Model setup
Two curvilinear mode grids were set up, west-grid and apex-grid (Fig. 1). The offshore boundary of both grids was set around the wave buoy B, of which the data were used to set up the offshore wave boundary. The west-grid encloses the entire North Sea coast of Sylt and the cross-shore extension is higher than 10 km. The highest grid resolution is about 15 m × 500 m (cross-shore × alongshore). The purpose of the west-grid was to transform offshore wave forcing up to the apex-grid and to avoid the wave shadow-zone effect at the lateral boundaries. The apex-grid covers the central part of the coast (~ 20 km), where there is the maximum curvature (i.e., the coastline orientation to the north and to the south from this domain hardly varies besides around inlets). The highest resolution is about 5 m × 450 m (cross-shore × alongshore). This grid was used to simulate the alongshore variability of storm erosion and sedimentation due to the coastline orientation.

Boundary forcing
From the observation period from January to April 2007, a simulation window was selected which consists of three severity events such that the first event has the highest storm severity, and the subsequent event has the lowest severity. This combination allows to investigate storm erosion under b. Grid setup a. Location c. Bathymetry different severity events and also the model skill in reproducing potential recovery after severe storm impacts.
The selected three events P1, P2 and P3 are shown in Fig. 2 for water level, wave height and wind characteristics. P1 shows the highest values of all these environmental forcings, which are synchronised, while they are the lowest in P2. The values in P3 are higher than P2, whereas significantly lower than P1. It should also be noted that both wave and wind approaches are from a westerly direction in P1 and P3. However, wave and wind directions appear to be opposite in P2. The storm power index was estimated using the approach of Dissanayake et al. (2015c), which defines by the square of storm wave height (Δh 2 ) multiply by duration (Δd) during the storm event.
Boundary forcing for D3D and XB within each storm period (P1: 7 days, P2 and P3: 10 days) was constructed using the data from the observations. Alongshore propagating tide was simulated by imposing the phase difference between the south and the north corner points of the offshore boundary. Astronomical tides at these two points were first estimated by the nesting approach in D3D between the calibrated German Bight model of Chu et al. (2013) and the setup domains. Finally, the water levels for the boundary forcing at the corner points were assessed by incorporating the estimated surge from the measured water levels at H (Fig. 1). Water level gradients along the lateral boundaries were set to zero in order to minimize the boundary disturbances. Offshore wave forcing was applied by spectral wave parameters (JONSWAP) from the buoy at B (Fig. 1). For the lateral wave boundaries in XB, the gradient of wave energy along the wave crest was set to zero to minimize the wave shadow-zone effect. A time-varying uniform wind field was implemented using the measured data from the land-based station at LW (Fig. 1). Characteristics of P1, P2 and P3 events are summarised in Table 2. Both environmental forcing parameters and the storm power index indicate that P1 is an extreme event, while P2 and P3 are weak and intermediate events respectively. Therefore, it can be expected highest storm erosion at the beach and dune system during P1 and the lowest impacts during P2.

Model simulations
First, a sensitivity analysis was carried out to investigate the effect of the lateral wave boundary location on the wave dynamics. D3D was simulated applying the apex domain alone, and the apex domain nested in the west domain. Both apex and west domains were separately simulated in XB. Wave dynamics from these five domains (3 in D3D and 2 in XB) were compared at four locations in the area of the apex domain.
Next, P1, P2 and P3 events were simulated in two series and the results were analysed using hydrodynamics and morphodynamics to identify the effect of coastline orientation on the storm impacts. In the first series, three events occurred in a cluster (i.e., the initial bathymetry for P2 was the storm impacted final bathymetry from P1 and so on). In the second series, events occurred in isolation (i.e., each event impacted on the same initial bathymetry).

Wave dynamics
Wave dynamics of the sensitivity runs are shown in Fig. 3 in comparison to the offshore wave characteristics (i.e., at B in Fig. 1). Four analysed locations are distributed in the apex-domain area. The first location (L1: 13 m depth) represents the domain middle, and the other three (~ 3 m depth) are along the shore (L2: south, L3: middle and L4: north). The offshore wave data indicate that the wave direction veers from West to North during the sensitivity runs, and the northwesterly wave has the peak storm wave height. According to the offshore wave direction, the wave shadow-zone effect is expected at the northern lateral boundary rather than at the southern. This is evident by comparing the D3D waves at L4 and L2, where the waves from the apex-domain alone (blue-line) show larger deviation at L4 than at L2 with respected the waves from the nested apex-domain (red-line). However, in the nested run, both apex-domain (red-line) and the west-domain (black-dash) show similar variations. L1 and L3, which are away from the later boundaries, have similar wave heights in all three domains. This further indicates the effect of lateral boundary on wave dynamics in D3D. In XB, both apex-(yellow-line) and west-domain (gray-dash) show similar variations of wave heights at all four locations implying that the notorious shadow-zone effect is hardly existing in XB. At the deep water point (L1), both models agree better, while at the shallow water points (L2, L3 and L3), waves in XB are higher than in D3D.
For the morphodynamic analysis, the apex-domain was therefore always simulated in D3D using the nested approach, while it was simulated alone in XB. Wave dynamics in P1, P2 and P3 were compared between models using the JONSWAP wave spectral density (Eq. 1), which is based on Hs and the mean period (Tm) (Hasselmann et al., 1973;Donelan et al., 1985;Nair and Kumar, 2017). XB predicts root-mean-square wave height (Hrms). Therefore, this was first converted into Hs following a relation of the Rayleigh distribution, Hs = √2×Hrms, which is however more applicable for normal wave conditions at deep water. The predicted wave data at the buoy L (see Fig. 1) were compared with the measured waves. Variations of normalised spectral density during the events are shown in Fig. 4 for both models with the measured wave heights at the B and L buoys. The spectral density was normalised with respect to the highest density in each event at B, so that P1, P2 and P3 can be compared together. For clarity, the spectral density is further shown for XB using Hrms. In P1, it appears that XB has higher density than D3D. However, the agreement between models is higher in P2 and P3 compared with P1.

Figure 4. Comparison of measured wave heights between B (red-) and L (black-line) buoys (Fig. 1) on the left panel, and the other panels: normalized JONSWAP wave spectral density using the predicted wave heights from Delft3D (D3D) and XBeach (XB: Hrms and XB: Hs) of the three events.
The agreement between the measured and the predicted wave data at the L buoy was analysed using the averaged wave spectrum over the events (Fig. 5). In the extreme event (P1), the measured data at B (13 m depth) and L (8 m depth) indicate that there is a one-peak spectrum towards low frequency (~ 0.1 Hz). This spectrum is fairly reproduced by both models, whereas the agreement between the data (L) and D3D is remarkable compared with that of XB. In the other two events, there are two-peak spectrums and the spectral density is one order of magnitude lower than P1. The main peak of the weak event (P2) still remains in the low frequency. This is captured by both models, particularly in XB (note, there is no second peak in D3D). The agreement between the models is higher than between the models and the data. In the intermediate event (P3), the main peak has been moved towards high frequency (~0.2 Hz). As found in P1, the spectrum of D3D better resembles the data compared with that of XB. The averaged spectral density curves were further compared between the data and the models using three statistical parameters; µ, σ and RMSE. µ is the relative difference, σ is the standard deviation and RMSE is the root-mean-square-error between the data and the models. The lower the parameter the higher the model skill in capturing the measured wave data. In all events, D3D claims the lowest value, particularly in P1 and P3 (Table 3). These results indicate that D3D has higher model skill in capturing nearshore wave heights than XB. Both models generally show, they can better produce high wave heights (P1, P3) compared with low wave heights (P2). The values of these statistical parameters further indicate, both D3D and XB can reasonably be used to investigate the nearshore wave dynamics under the default parameter setting.

Morphodynamics
The alongshore variation of storm impacts under different storm events and combinations was analysed using cross-shore profile evolution, bed level and sediment volume change, and alongshore erosion and sedimentation pattern. The analysis area for morphodynamics covers the beach and dune system from -5 m to +5 m MSL of the apex-domain (see red-enclosed area in Fig. 1).
Cross-shore profile evolution during the impacts of storm events in Cluster (I) and Isolation (II) is shown in Fig. 6 for the six profile locations from the South to the North (C1 to C6) of the central coast. Both the storm applications and the models caused strong erosion at the beach and dune system, and then sedimentation in the nearshore area. P1 has resulted in the highest impacts along the coast. In Cluster, D3D shows lower erosion (particularly above MSL) than XB. The difference among the storm events in D3D is marginal at the southern and the northern profiles. However, the profiles around the apex of this coast (C3 and C4) show increase of erosion (close to MSL) from the subsequent events. In XB, the dune erosion is prominent at all profiles. After the first extreme event, the other events caused further erosion around MSL along the coast. Both models show relatively higher erosion and sedimentation around the apex than at the South and the North of the coast. In Isolation, erosion and sedimentation occurred based on the severity of each event. Dune erosion (above +2 m MSL) is found only during the first event. Thereafter, the profile change is marginal in both models along the coast. The area around the apex still shows high morphodynamics from both models, of which XB resulted in higher dune erosion than D3D. Bed level change (upper panels in Fig. 7) indicates that XB resulted in strong erosion above the dune toe level (DT). After the first extreme event, there are hardly any changes of bed levels in P2 and P3 in Cluster. However, the latter two events caused erosion below mean high water level (MHW) during the impacts of events in Isolation. The majority of bed level change in P2 and P3 is within the 20% limit compared with the no bed level change situation. Fig. 1

) in each event from Delft3D (D3D) and XBeach (XB) during the application of Cluster (I) and Isolation (II). Dash-line of bed level change indicates 20% deviation from no-change situation.
Sediment volume change shows the amount of erosion (-) and sedimentation (+) in the different depth classes within the analysis area (lower panels in Fig. 7). The first depth class (-5 to -3 m MSL) experienced the highest erosion and the highest sedimentation in both models. The sedimentation amount is twice of the erosion because of the eroded sediment approaching from higher elevations of the beach and dune area. From MSL towards dune, the erosion from D3D decreases while increasing in XB (see 3 to 5 m depth class). As found with the bed level change, the amount of erosion and sedimentation of the individual events is higher in Isolation than in Cluster. In both cases, XB caused higher amount of erosion than D3D.
Alongshore variation of the normalised erosion and sedimentation pattern is shown in Fig. 8 for P1, P2 and P3 in Cluster and Isolation. Results are shown for the coastal stretch within the analysis area (see Fig. 1). Values were normalised by dividing with the highest erosion volume of the P1 event. With respect to the reference line (noted by 0), landward shadow zones indicate erosion and seaward is sedimentation. Net alongshore variability is shown by the dash lines with the respective colour for each model (D3D: blue and XB: red). In the extreme event (P1), XB produced the largest erosion and sedimentation along the coast. The D3D results are generally low, whereas they follow the alongshore variability of XB. Around the apex, both models resulted in higher erosion than the other parts of the coast. The weak event (P2) has almost no impact on the coast in Cluster and minor changes in Isolation. The latter still indicates that the area around the apex has experienced strong erosion and sedimentation. In the intermediate event (P3), both models show almost similar morphodynamics in Cluster and Isolation from the apex to the north of the coast (note. there is no dune erosion). The net change of alongshore variation shows that strong erosion at the apex in all events, but different amounts.

DISCUSSION
D3D and XB have different approaches to evaluate the propagation of offshore wave spectrum (JONSWAP). D3D uses SWAN, which estimates the evolution of offshore short-wave spectrum in both frequency and directional spaces (Booij et al., 1999), while XB computes the long-wave envelop in the short-wave group scale using a mean frequency in the directional space only (Roelvink et al., 2018). In XB, the gradient along the wave crests of the wave energy can be set to zero. This allows the crests of the wave groups to have the correct orientation along the lateral boundaries. In D3D-SWAN, there is no wave energy along the lateral boundaries, unless the user defines waves. Therefore, the shadow-zone effect is inevitable in D3D, which needs to be avoided applying large domains or model-nesting. However, D3D-SWAN has high skill in wave prediction being a fully spectral (frequency and direction) wave model. This was found with the wave dynamics at the L buoy, waves from D3D agreed better with the data than from XB. It should be noted that the predicted Hrms of XB was converted to Hs using a relation (Hs=√2Hrms) based on the Rayleigh distribution, which is however more applicable for normal wave conditions at deep water. Thus, the agreement of Hs with the data was higher in P2 (weak waves) than in P1 (strong waves) in XB. During storms, short-wave dissipates in the surf zone and long-wave approaches the coast. For the progressive dune erosion, long-wave impact, avalanching and return flow are important processes. These physical processes are embedded in XB (see Table 1). Storm erosion showed that XB resulted in higher erosion along the coast than D3D. Therefore, XB is suitable to estimate the beach and dune erosion during extreme events. However, during calm events, short waves dominate, which indicates the relevance of D3D in simulating calm events (Van Ormondt et al., 2020). This is evident with the results by the fact that during the weak event, the difference of erosion between models is low (see Fig. 6 and 7).
Default parameter setting of both models resulted in similar trend of storm impacts along the coast. Wave dynamics showed that there are reasonable agreements with the measured data (µ: -0.01-0.63, σ: 0.23 -0.68, RMSE: 0.13 -0.45 m 2 /Hz). Storm erosion in the models had the similar pattern of alongshore variation. However, the occurrence of the weak event (P2) after the extreme event (P1) could not produce any recovery at the beach and dune system. This indicates the requirement of different parameter setting to promote the onshore transport. The user has to adjust the parameters mainly related to wave asymmetry and skewness in both models to decrease offshore-and to increase onshore-transport (e.g., Pender and Karunarathna, 2013).
Storm erosion along the coast showed the sensitivity to storm severity, nearshore morphology and the storm landfall, which depends on the coastline orientation and the storm approach. Selected three events had different severity levels and they resulted in distinct erosion and sedimentation at individual locations along the coast, which were proportional to the storm power index. After the first extreme event, the nearshore morphology of the subsequent events was weaker in Cluster than in Isolation. Therefore, the morphodynamic impacts of the last two events were strong in Isolation. Similar findings were discussed in Dissanayake et al. (2015c). The first event (P1) approached from the SW -NW sector (wave from NW and wind from W), and was severe to impact the entire coast (see Fig. 6 and 8). The convex-shape of the coastline provided however different storm landfall at the individual locations along the coast. This is evident by the storm erosion, which showed higher impacts along the stretch from the apex to the North than at the southern coast. Furthermore, the highest erosion of this event was found at the apex (see net change in Fig. 8). In the last event (P3), both wave and wind approached from NW. Resulting erosion pattern was again strong from the apex to the northern coast, while showing the highest values at the apex. Therefore, high storm impacts and erosion are expected from the apex to the North, and the apex has the highest vulnerability of this convex-shape coast, when storm events approach from the northwesterly directions. The dominant wave direction of this coast is from NW. Thus, the latter coastal stretch draws attention in terms of the application of management strategies (e.g., large quantities of sand nourishment than the southern coast). Furthermore, on the convex-shape Sefton coast (Liverpool bay, UK), extreme events (wave and surge) occur from SW to NW winds (Brown et al., 2010). This coast has also experienced strong erosion around the apex (Esteves et al., 2011). Therefore, it appears, the beach and dune system at the apex is more susceptible to storm impacted morphodynamics on the convex-shape coastlines.

CONCLUSIONS
A hybrid numerical modelling approach using Delft3D and XBeach was adopted to investigate the effect of coastline orientation on storm erosion on the convex-shape coastline, the North Sea coast of Sylt. Both models were simulated (with default parameter setting) applying three different severity events (extreme, weak and intermediate in sequence) considering two occurrences of events (cluster and isolation). The model predicted wave dynamics and morphodynamics were analysed to evaluate the storm impacts.
Statistical comparison of wave dynamics showed that there are reasonable agreements with the measured data (e.g., RMSE of the wave spectral density ≤ 0.45 m 2 /Hz). Waves of Delft3D better resembled the data. In the weak event, the difference of wave spectrums (JONSWAP) between the models is low. The conversion from the XBeach predicted waves (Hrms) to Hs is generally valid for calm waves. Severity of the events and the initial morphology determined the storm erosion. The extreme event had the highest erosion, while the events in isolation resulted in higher erosion than in cluster. Alongshore variation of storm erosion depended on the coastline orientation and the storm approach. In the extreme event, wave from NW and wind from W approached the coast. The convex-shape of the coastline provided however different storm landfall at individual locations along the coast. The coastal stretch from the apex to the North showed higher erosion than at the southern coast. In the intermediate event, both wave and wind approached from NW and produced high impacts along the same stretch as in the extreme event. In both events, the highest erosion occurred around the apex. In the last two events, the alongshore variation of storm erosion was fairly similar. With respect to the model numerics and the predicted morphodynamics, XBeach is suitable to investigate the impacts (beach and dune) during extreme events, while Delft3D is appropriate during normal weather conditions. This coast endures storm events from the northwesterly directions. Therefore, the apex (particularly) and the northern coast are more vulnerable to storm impacts than the southern coast. These findings agree with the present management strategies, which appear that sand nourishments on this coastal stretch are larger than the southern coast to withstand the storm erosion.
These results however provide the first insights of the effect of coastline orientation on storm erosion along the convex-shape North Sea coast of Sylt. Further studies are prerequisite in order to quantify the alongshore erosion with respect to the measured event-bounded morphodynamics.