The impact of wind‐waves and sea level rise on the morphodynamics of a sandy estuarine shoal

Intertidal shoals are pronounced morphological features found in many estuaries worldwide. Apart from maintaining an ecologically unique intertidal environment, shoals also protect adjacent dyke systems by attenuating waves. The fate of sandy shoals under anticipated sea level rise (SLR) scenarios is underexplored.

boundaries, while other shoals exist in a free form surrounded by deep tidal channels (Hibma et al., 2004;Robinson, 1960).
In a tide-dominated, sandy estuary, typical morphological features are ebb and flood channels meandering between shoals (Ahnert, 1960;Robinson, 1960). The channel-shoal patterns arise from positive feedback between tidal currents, sediment transport, and initial bedforms (Hibma et al., 2004;Schramkowski et al., 2002;Schuttelaars & Swart, 1999;Seminara & Tubino, 2001), while the plan form of the estuary has an impact on their location and bar scales (Dam et al., 2016;Lanzoni, 2002;Van der Wegen & Roelvink, 2008, 2012. The intertidal morphology is governed by a complex interaction between tidal hydrodynamics, wave action, and sediment supply (e.g., Friedrichs, 2012;Kirby, 2000;Le Hir et al., 2000;Roberts et al., 2000). Tidal forcing usually causes onshore sediment transport, while waves enhance shoal erosion and prevent deposition (Christie et al., 1999;Le Hir et al., 2000;Yang et al., 2003). In response to increased sediment supply and a large tidal range, convex-up cross-section profiles are likely to form, while concave-down profiles are associated with decreased sediment supply and increased wave action (Bearman et al., 2010;Friedrichs, 2012;Hunt et al., 2015). Wind-waves play an important role in morphodynamic states (Fagherazzi et al., 2007;Maan et al., 2018;Van der Wegen et al., 2019). High wave energy tends to generate a steeper intertidal profile Friedrichs, 2012). By studying the Dutch Wadden Sea and Eastern Scheldt, respectively, Janssen-Stelder (2000) and Kohsiek et al. (1988) highlighted that a combination of tidal currents and wave action is responsible for the erosion of the intertidal flat during storm conditions. DeVet et al. (2018) showed that a sandy shoal may propagate along the wave direction and wind direction. Elmilady et al. (2020) highlight that wave-induced resuspension is essential for distributing sediment supplied from deep channels over the intertidal flat.
Sea level rise (SLR) will potentially have considerable impacts on estuary morphodynamics. Recent studies applying morphodynamic numerical models for estuarine systems show that sea level rise can lead to noticeable intertidal area loss (e.g., Dissanayake et al., 2012;Elmilady et al., 2019;Mariotti & Fagherazzi, 2010;Van Der Wegen, 2013;  and a steeper slope of tidal flats (Mariotti & Fagherazzi, 2010). SLR-induced increased water depth weakens wave action on channel bank, where sediment deposition is favored. Despite sedimentation, intertidal flats cannot keep pace with an exponentially rising sea level, leading to a loss of intertidal area along with increasing inundation frequency (e.g. Elmilady et al., 2019;Elmilady et al., 2020;Ganju & Schoellhamer, 2010;Van Der Wegen, 2013; . Van Goor et al., (2003), Lodder et al. (2019), and Elmilady et al. (2020) suggest that intertidal flats can be maintained under a limited, linear SLR after a decadal adaptation time during which intertidal area is lost. This maintenance only exists when the accretion rate of the intertidal flats reaches the SLR rate. This usually requires a long adaption timescale and depends strongly on the local conditions of tides, waves, and sediment supply (Elmilady et al., 2020;Van Goor et al., 2003).
This synopsis reveals some intriguing knowledge gaps. Tidal basin models covering a large domain often use a coarse resolution of above 100 m to save computation time (e.g., Dam et al., 2016;Hibma et al., 2004;Van der Wegen & Roelvink, 2008;Wang et al., 2007), whereas smaller tidal flat models have a much higher resolution of typically 10-50 m (e.g., DeVet et al., 2018;Maan et al., 2018). The coarse resolution could disregard the subgrid (wave) dynamics and morphological structure, while the higher resolution is time-consuming. In addition, while tidal flat models stress the importance of wind-waves, many tidal basin models do not include wave action to reduce modeling complexity and computation time (e.g., Dam et al., 2016;Van der Wegen & Roelvink, 2012;Wang et al., 2007) or include waves which only impact mud transport (Braat et al., 2017). Tidal flat models including waves focus on the mud, while many estuarine shoals are sandy (e.g., Western Scheldt (Robinson, 1960), Waddenzee (Wang et al., 2012), and the Dyfi Estuary (Brown & Davies, 2010)). Furthermore, many studies apply a schematized and simplified model setup (e.g., Braat et al., 2017;Van der Wegen & Roelvink, 2012), disregarding site-specific conditions such as dredging and disposal activities for navigational purposes. Finally, although sea level rise threatens the survival of estuarine shoals and flats, limited processbased morphodynamic modeling studies have considered the impact of wind-waves when investigating anticipated sea level rise scenarios.

| AIM AND APPROACH
This research aims to explore the impact of wind-waves on the longterm morphodynamic evolution of estuarine sandy shoals in a realistic setting including the impact of sea level rise.
We apply a 2D process-based numerical model (Delft3D) to simulate the evolution of a shoal complex located in the landward side of the Western Scheldt estuary (the Netherlands). First, we assess the skill of our model in simulating decadal morphological development in a complex estuarine environment by qualitatively comparing modeled and measured bathymetries. We then compare the modeled shoal evolution over 50 years with and without waves and investigate the underlying reasons for the morphological difference caused by waves. Finally, starting from the 2013 modeled bathymetry and keeping the same wave environment in each scenario, we perform a 100-year forecast by imposing different sea level rise rates and study the impact of wind-waves on shoal evolution under SLR.

| STUDY AREA
The Western Scheldt estuary (the Netherlands; Figure 1) is one of the youngest natural estuaries in Western Europe. The large estuaries and tidal basins in the Netherlands were formed around AD 200 (Bosboom & Stive, 2012). Storms, floods, subsidence, and a continuing small sea level rise led to a rapidly increasing size of the Western Scheldt in the late Middle Ages ( Van der Spek, 1994). Over time, various smaller secondary embayments have been reclaimed, decreasing the area of the Western Scheldt and with some impacts on the morphology of the main channels (Dam et al., 2016;Nnafie et al., 2019). Nowadays, it is a funnel-shaped mouth of the River Scheldt. The tide propagates up to Ghent, which is located around 160 km upstream (Figure 1). Its cross-section width decreases roughly exponentially from the estuary mouth at Vlissingen (around 5 km wide) to the estuary head near Ghent (less than 100 m wide). The width-averaged depth decreases from about 15 m at Vlissingen to 3 m at Ghent. The morphology of the Western Scheldt is a multiple-channel system with ebb channels and flood channels separated by intertidal shoals (Van Veen, 1950). It is a sand-dominated estuary with less than 10% mud, mostly found alongside the estuarine margins and salt marshes (Kuijper et al., 2004).
Major dredging activities have taken place since the 1970s to guarantee navigation to the Port of Antwerp (Vroom & Schrijvershof, 2015;Wang et al., 2002).
The estuary is a tide-dominant system that is subject to strong spring-neap tidal variations. The most dominant tidal components are semi-diurnal (M2) and quarterly diurnal tides (M4) (Wang et al., 2002).
The mean tidal range is 3.8 m at Vlissingen, 5.2 m at Niel, and 1.9 m at Ghent, with ratios between the rising and falling durations of 0.88, 0.75, and 0.39, respectively. The annual-averaged discharge from the River Scheldt at the confluence of the Rupel and Scheldt amounts to about 110 m 3 /s. The yearly variation ranges from 50 to 200 m 3 /s (Kuijper et al., 2004). In the landward portions of the Western Scheldt, local winds generate waves with an average significant wave height of 0.2 m, propagating mainly from the southwest, northeast, and northwest.
This study focuses on the shoal complex (right panel in Figure 1) including the Shoal van Ossenisse, the Molenplaat, and the Rug van Baarland located in the mid-estuary, about 30 km landward from the mouth. This area is chosen because it is subject to limited penetration by North Sea offshore waves, which allows for the implementation of a wave model covering a part of the inner estuary and only considering small locally generated wind-waves. Shoals and the adjacent channels in this area are sand dominated with D 50 ranging from 120 μm to 250 μm (Kuijper et al., 2004).

| MODEL CONFIGURATION
In this research, we apply a two-dimensional horizontal (2DH) Delft3D model to simulate the long-term morphological evolution of the Western Scheldt estuary. This is a coupled flow and wave model ( Figure 1). The flow model solves the shallow water equation and includes morphodynamic changes. The sediment transport formula distinguishes suspended load transport and bedload transport according to a reference height (Van Rijn, 1993). The wind-generated waves are calculated by the spectral wind-wave model SWAN (Booij et al., 1999). This describes waves by the two-dimensional wave action density spectrum. The flow model is coupled to the wave model through online coupling and the calculation of the currentwave interaction. More descriptions of the model terminology are reported in the Supporting Information.
The implemented model configuration is based on the NeVla-Delft3D model, which is calibrated for hydrodynamic conditions (Maximova et al., 2009; and has been previously  Sutherland et al. (2075), the BSS is calculated to quantify the model skill using Equation 1. The skill score shows that the landward part performs better than the wider mouth and about 45% of the modeled erosion and sedimentation volume in the landward domain meets the criterion BSS ≥ 0.5.
in which Δvol is the volumetric change compared with the initial bed (m 3 ), mod is the modeled quantity, meas is the measured quantity, and brackets denote the arithmetic mean.
We divide the NeVla flow model grid into five subdomains to allow for parallel processing and local grid refinement. The seaward boundary extends 20 km from the estuarine mouth while the land-  Table 1. This modeling configuration takes about 3 days to simulate 6 months of hydrodynamic time on an eight-core (2.6 GHz) machine.  The aim of the current work is not to reproduce observed developments but rather to explore the impact of wind-waves on morphodynamic development. We consider the model accuracy to be high enough and to generate a stable morphology in the surroundings of Shoal van Ossenisse to continue exploring the impact of windwaves on morphodynamic development.

| Impact of wind-waves on shoal evolution
We investigate the impact of wind-induced waves on the intertidal area by comparing simulations with waves (Wave) and without waves (NoWave). Figure 4 shows the bathymetry difference (Wave -NoWave) over 50 years when applying different wind directions.
Red/blue color indicates a higher and lower bed level for the Wave than NoWave case. Note that higher/lower bed level can be due to more/less sedimentation or less/more erosion for the Wave than  an estuary (Boon, 1975). Figures 6a, 6c,

| Processes underlying wind-wave impact
In this section, we substantiate an explanation for the elevation lowering and widening of the intertidal area along with the channel widening and deepening at large depths. We hypothesize that waves erode the intertidal area by enhancing the bed shear stress on the shoal. This raises suspended sediment concentrations in the water column at the shoal. Combined with wave asymmetry and wind-driven flow, a high SSC enhances sediment transport from the shoal towards the channel in the direction of wave propagation and causes higher prevailing SSC on the leeside of the shoal. In addition, increased SSC levels will occur at the shoal edges by ebb currents draining the shoal. As a result, the shoal widens and the channel narrows. This leads to flow convergence and erosion in the channel with higher velocities, deepening the channel.
To investigate the sediment transport variation due to waves, we ran scenarios on a fixed bathymetry (without morphodynamic updates). Model results show that intratidal variations (from Supporting Information, Figure 2h) of the bed shear stress and SSC have peaks when the water level is around MSL. At this stage, the large flow velocity enhances nonlinear wave-current interactions, increasing bed shear stress and inducing high SSC, especially on the upwind side. During low water slack, the influence of waves can be extended to the subtidal domain, but the amount of area feeling wave motion is limited. During high water slack, the water depth is too large, so that the bed level hardly feels wave motion anymore. Tidal residual bedload sediment transport is more directed towards the wave propagation direction (Figure 9c). This is because, in Delft3D, the bed load transport equation includes the suspended sediment transport caused by wave asymmetry, which represents a large portion of the total bedload sediment transport over the shoal.
However, the bedload sediment transport magnitude is negligible compared with the suspended sediment transport magnitude and barely contributes to the total sediment transport. NoWave and Wave in both NoSLR (Figure 11a and 11b) and SLR (Figure 11c and 4d) shows that the wave impact on shoal patterns becomes more pronounced in the following 100 years. In particular, the prescribed single wave climate causes continuous shoal migration in the downwind direction ( Figure 11g). The SLR scenario sees faster propagation due to wave actions (cf. Figure 11g and 11h, which show more sediment at the shoal leeside in Wave when imposing SLR 1.10 m). This is consistent with Elmilady et al. (2019), who showed that, as SLR increases the water depth on the shoal, the shallower area of the shoal starts to experience greater wind-wave influence.
This increases the landward extent of wave attack on the upwind side and therefore sediment on the downwind side.
The study area experiences net erosion with a slightly increased channel volume and a decreased intertidal area and intertidal volume over 100 years ( Figure 12). Figure 12e shows that SLR enhances the decrease of the total sediment volume, but the decrease is not a linear function of SLR rate and the trend line is not smooth but rather irregular and sometimes shockwise. For example, in NoWave cases, linear and nonlinear 1.10 m SLR, linear and nonlinear 1.95 m SLR lead to 9%, 47%, 62%, and 91% more erosion than NoSLR with respect to the initial state, respectively. Figure 12a shows an overall decrease of the intertidal area in the estuary. The development of the intertidal area does not keep pace with the rising sea level due to the SLR-induced increased inundation. The higher the sea level rise, the greater the decrease of the intertidal area. The nonlinear and linear SLR signals show a small difference. Initially, the decrease of intertidal area for nonlinear SLR signal is milder. Afterwards, the decrease rate becomes larger due to the faster inundation rate. Similar to the intertidal area, SLR causes a notable drop in the intertidal volume (Figure 12b), while the decrease does not linearly relate to SLR rate. A significantly larger channel volume (Figure 12c) by about 13% and 23% compared with NoSLR arises from 1.10 and 1.95 m SLR, respectively, by 2100. This is partly because of the increase in the inundated area below MLLW.
Additionally, the channel volume calculated below the constant initial MSL (Figure 12d) suggests that SLR causes net erosion of the channel with a percentage of around 3% and 7% larger than NoSLR for the signal with 1.10 and 1.95 m.
The impact of waves does not differ much among the different SLR scenarios. Compared with the SLR impact, the difference in channel volume due to waves (less than 1%) remains limited. This implies that the magnitude of wave impact is much less than the magnitude Wave action on modeled channel-shoal patterns shows some recognized morphological impact. We observe that shoals tend to migrate in the direction of wave propagation, which is consistent with results by DeVet et al. (2018). Furthermore, the larger intertidal area resulting from waves hence leads to channel deepening because of velocity convergence.
Our approach implicitly assumes that developments and interactions with the surroundings of the shoal (i.e., the larger system) are of secondary importance to the processes that we wanted to study (wind-wave impact and SLR action. Therefore, the estuary is seen to react more to SLR than waves, and including waves or not does not change the estuarine reaction to SLR. This suggests that, if computational cost is a limitation, the process-based morphodynamic modeling could focus more on the SLR impact than small wind-waves.

| CONCLUSIONS
The current research explored the long-term morphodynamic evolution of estuarine sandy shoals under the impact of small wind-waves and sea level rise by means of a numerical, process-based model (Delft3D) in a realistic setting. Firstly, we run models for 50 years, only including wave as a variable. Later, we imposed SLR at the seaward boundary and run models with and without waves under different SLR scenarios for 100 years.
Over the decadal timescale, we observe morphological differences caused by small wind-waves. The estuarine shoals tend to migrate slightly along the direction of wave propagation. The top of the shoal is lowered, with more sediment at shoal edges and downwind side. The channel is wider at larger depths and deeper in the bottom, especially near the areas of wider shoal edges caused by waves.
We investigate the underlying mechanism behind these observed differences. Waves erode the intertidal area by enhancing the bed shear stress on the shoal. The sediment resuspended by waves is partly transported to the leeside of the intertidal area and partly occurs at the shoal edges by ebb currents draining the shoal. As a result, the shoal widens with a decreased channel area. This leads to flow convergence, and deep parts of the channel widen and deepen in response.
However, waves do not lead to a fundamental difference in estuarine evolution. On a long-term timescale, waves keep causing shoal migration in downwind direction and lower the shoal, while the quantitative differences of intertidal area and volume, and channel area and volume diminish or stay constant. On the long-term estuarine scale, small wind-waves are less dominant than the other governed factors that are mainly responsible for the readjustment of channelshoal patterns, such as tidal asymmetry. The presence of wind-waves does not fundamentally change the general evolution of the estuarine system. The impact of SLR dominates the effect of wind-waves on the estuarine morphological evolution. The shoal elevation keeps growing with the rising sea level rise, while a notable intertidal area and volume loss still occurs due to the increased inundation. The combination of the increased submerged channel volume and a net channel erosion leads to a larger channel volume caused by SLR. Waves retain the function of lowering the intertidal area and causing channel profile variations, but their impact is limited compared with SLR and the general dynamics of the estuary. The inclusion of waves does not fundamentally change the reaction of the estuary to SLR.
We have different proposals for whether or not to include small wind-waves in the process-based numerical model when studying the estuarine morphology. If one investigates the development of the intertidal flats, we recommend the inclusion of highresolution small wind-waves on shoals, as wind-waves are major drivers of the intertidal morphodynamics, especially in environments with high wave energy and fine sediment. However, lowresolution models without wave actions could be used to study long-term and large morphodynamic development of a sandy, tidedominant estuary, if the grid resolution adequately captures the tide-residual sediment transport patterns. The minor impact of wind-waves compared with SLR on the estuarine morphological development on long-term timescales suggests the insignificance of including small wind-waves when investigating anticipated sea level rise scenarios.
Process-based models include many interacting processes and encompass multiple timescales, which makes it difficult to discriminate various impacts. In addition, the extensive dredging activities in the channel and continuous deposition of dredged materials on the shoal may influence the channel profile and disturb the analysis.
This work reveals the significant impact of wind-waves on the morphology of sandy shoals, consistent with literature that studies muddy intertidal flats. Nevertheless, our findings show the minor impact of small wind-waves on a sandy estuarine system as a whole, which is in contrast with some mud studies. Since muddy systems are more sensitive to waves and have different sediment supply sources, they may show distinct behavior. In reality, most estuaries are mixed with sand and mud, thus future research should focus on implementing our approach to a system with multiple sediment fractions.