Vegetation Reconfigures Barrier Coasts and Affects Tidal Basin Infilling Under Sea Level Rise

Worldwide, many tidal basins associated with barrier coasts have infilled over the past millennia due to the combination of sediment supply, wave‐tidal sediment transport, and eco‐engineering effects of vegetation. However, the biogeomorphological interactions between saltmarsh and the morphodynamics of an entire coastal barrier system are poorly understood, especially under sea level rise (SLR). Here, we study the evolution of a barrier coast for combinations of mud availability, presence of vegetation, and SLR. We developed a novel biogeomorphological model of an idealized barrier coast enclosing a tidal basin with sandy‐clayey sediments that was subjected to tides and waves for a century. The morphodynamic Delft3D model was coupled to a vegetation code which accounts for the dynamics of marsh‐type vegetation. Initially, vegetation contributed to reducing the tidal prism while sediment was imported. However, with SLR this trend was reversed and the tidal basins started to export sediment for vegetated runs after about 50–60 years while the unvegetated scenarios continued to infill in pace with the SLR. The sediment export was caused by cascading biomorphodynamic feedback effects triggered by vegetation which modified channel and shoal dynamics. Even under higher mud supply, the SLR resulted in vegetation collapse. The hypsometries, similar to natural systems, showed that vegetated systems converge to an alternative stable state condition. We conclude that the long‐term resilience of the tidal basin associated with sediment infilling under SLR can be reduced by cascading large‐scale effects of vegetation on the morphodynamics of barrier coasts.


of 22
combined with sediment transport processes exert a major control on the formation of channels, the ebb-flood delta and on the inlet dimension and thus shape the entire barrier coastal system (De Swart & Zimmerman, 2009).
Many tidal basins worldwide were formed during sea level rise in the early and middle Holocene (Boyd et al., 1992;de Haas et al., 2018;FitzGerald & Buynevich, 2003). Tidal basins receive sediments from a combination of coastal and fluvial supplies. While a high fluvial supply can be present for estuaries, the fluvial supply in tidal basins can be absent or negligible when the open coast is the main source of sediments. The coastal sediment supply derives from the advection of mud, littoral drift of sand, and ebb delta dynamics. Under increasing sediment demand, for example, due to sea level rise, the basin needs more sediments to keep up with the increase in accommodation space. Sediments can be supplied by the ebb delta in the short term while in the long term the offshore and the updrift coast, including the barrier itself, are the main sources (Beets & van der Spek, 2000). Mud and sand show different behaviours. Noncohesive (coarser) sandy sediments are transported as bed and suspended load and they respond strongly to the flow velocities due to the nonlinear relation of flow velocity and sediment transport (van Rijn, 2007a). Cohesive muddy (finer) fractions are transported as suspended load and their deposition is strongly related to the duration of rising-falling tides and the slack-water periods that allow for sedimentation. Sandy material is supplied by tidal and wave-related sediment transport from the adjacent barriers and ebb-tidal delta. Mud and fine suspended material advected from offshore and stirred from shallow areas within the basin are transported by wave-tidal currents. Part of these sediments are bypassed to the downdrift coast via the inlet, ebb delta, and from offshore, while another part remains inside the basin (Bruun & Gerritsen, 1959;Elias et al., 2019;FitzGerald, 1996;FitzGerald et al., 2000;Lenstra et al., 2019). The balance between imported and exported sediments is influenced by the tidal asymmetry, the shape and length of the basin, intertidal volume storage, the antecedent hypsometry (bed elevation distribution of channels and shoals), and the presence of vegetation (Friedrichs & Aubrey, 1988;Stark et al., 2017;van Maanen et al., 2013).
The net sedimentation (infilling) is therefore a result of the hydrodynamics, sediment transport, and eco-engineering effects of vegetation. Eco-engineering means the capacity of vegetation to alter flow and sediment transport processes on biomorphodynamics. Paleogeographical reconstructions of the Netherlands show that middle Holocene back-barrier systems were filled in by reed peat (Vos, 2015) and transformed into a coastal plain when the sea level rise rate decelerated from ca. 5 mm/yr to less than 2 mm/yr (Hijma & Cohen, 2019). This suggests a long term sustained feedback of vegetation promoting sedimentation and consequently tidal basin filling . Modelling studies mostly support that vegetation retains more sediment and contributes to the local increase of bed levels (e.g., D'Alpaos, 2011;Fagherazzi et al., 2012;Kirwan & Murray, 2007) and that for sea level rise exceeding the sedimentation due to sediment supply, the vegetation can collapse, exacerbating the drowning (Kirwan et al., 2010;Schuerch et al., 2018). Other studies suggest that large patches of vegetation reduce flow velocities, local sediment supply, and sedimentation away from the feeder channels (Boechat Albernaz et al., 2020;Brückner et al., 2020;D'Alpaos et al., 2007).
Mud and vegetation tend to confine channels (i.e., smaller width to depth ratio) which favors ebb-dominance (Bij de Vaate et al., 2020;Friedrichs & Aubrey, 1988). Depending on the local small-scale characteristics, adjacent channels and shoals can be flood or ebb dominant, whereas on average the channels are ebb dominant and shoals flood dominant (van Veen et al., 2005) due to the balance between friction and inertia. Shallow shoals and intertidal flats also allow for the settling of vegetation (assuming that the underlying environmental conditions are also satisfied), which in turn alters the hydrodynamics and sediment transport in the morphodynamic feedback loop (Boyd et al., 1992;Wright & Thom, 1977).
Tidal basins depend on sediment fluxes to/from the littoral zone, and the evolution of each subsystem affects the interaction between the parts (Robbins et al., 2022). The resulting sediment transport in and out from the basin depends on the wave-tidal dynamics that, in turn, depend on the tidal basin and ebb-flood delta properties, including effects of vegetation. The working hypothesis is that the tidal basin development is the result of interactions between the littoral sandy processes, the mud dynamics and vegetation. As such, the sediment fluxes cannot simply be assumed or imposed at the basin boundary (e.g., inlet) but are the dynamic result of coastal processes interacting with the tidal basin.
The combined effects of vegetation changing local sedimentation, channel configuration, tidal penetration, and distal sedimentation suggest a need to evaluate the conditions under which vegetation promotes a negative feedback (increasing infilling towards a steady state) or a positive feedback (erosion and drowning) on the basin evolution, regardless of the sediment supply. These coupled mechanisms have not yet been investigated with 10.1029/2022JF006703 3 of 22 more comprehensive hydro-biomorphodynamic models, while the negative feedback was inferred from paleogeographical reconstructions .
In such a complex open biogeomorphic system, it is an open question under which conditions a steady state, or equilibrium, occurs. Research on barrier coasts and tidal basins started with field-derived equilibrium conditions attributed, for example, to inlet dimensions and tidal prism (Escoffier, 1940(Escoffier, , 1977Huismans et al., 2022;O'Brien, 1967;Powell et al., 2006) or the volume of ebb and flood delta related to basin dimensions (Dronkers, 1998;Powell et al., 2006;Walton & Adams, 1976). This data-driven approach is valuable when studying observation-rich systems that are in a steady state. However, they provide less insight about the underlying physical processes, especially under changing boundary conditions (Ranasinghe, 2020). Numerical and analytical models have been applied to overcome these limitations as they can test alternative scenarios including changes in hydrodynamics and sediment supply. However, comprehensive process-based models, such as Delft3D (Lesser et al., 2004) and TELEMAC (Villaret et al., 2013) are complex and computationally intensive, especially with wave-driven sediment transport. For these reasons, morphodynamic models until now were able to perform either short timescale simulations (O-years) or needed to be simplified (e.g., no waves, single sediment fraction or neglect the coast and inlet) to simulate long-term morphodynamics (O-decades). Only lower complexity models were able to reach the time-space scales of full barrier coast development (De Vriend, 1991;Elias et al., 2019;Ranasinghe, 2020). Furthermore, the effects of eco-engineering vegetation species  and eco-engineering benthic species (Brückner et al., 2021), often disregarded in models, have large eco-morphodynamics effects. In summary, we currently rely for our long-term predictions on a priori equilibrium modes or simplified models (D'Alpaos, 2011;D'Alpaos et al., 2007;De Swart & Zimmerman, 2009;Fagherazzi et al., 2012;Kirwan & Murray, 2007;Leuven et al., 2019;Mariotti & Canestrelli, 2017;Nardin et al., 2020;Ranasinghe, 2020;Townend et al., 2016;Xu et al., 2019).
In view of the raised questions and ongoing global climate change, there is an urgent need for understanding how barrier coasts will respond to natural and human induced relative sea level rise and changes in sediment supply (e.g., Dunn et al., 2019;Eslami et al., 2019;Nienhuis et al., 2020;Syvitski et al., 2005). The starting point in this paper is the hypothesis that the interaction between sediment supply and vegetation determines the evolution of barrier coasts under rising sea level. We expect that the local morphodynamic feedback between sediment transport and vegetation on the formation of channels and shoals exerts an important control on the basin-scale evolution.
In order to test our hypotheses and unravel the main underlying mechanisms behind barrier coast evolution and basin infilling, we developed a novel large-scale biomorphodynamic model consisting of an entire idealized barrier coast system with a dynamic vegetation module that interacts with the hydro-morphodynamics. The domain encompasses a tidal basin with a single inlet within an alongshore uniform coast. The model is the first process-based (i.e., hydrodynamics resolving model coupled to detailed sediment transport formulations) biomorphodynamic model that includes comprehensive wave climate, morphological tidal constituents including overtides, a sediment composition from sand to clay, and dynamic vegetation. Our scenarios cover a century of biomorphodynamics under sea level rise for a range of concentrations of offshore mud supply. In total, we simulated eight scenarios with and without vegetation, with different clay concentrations between 15 and 60 mg/L and reference scenarios without sea level rise.

Methods
The morphodynamic simulations were performed in a depth-averaged (2DH) configuration of Delft3D (Deltares, 2020, modified from tag 7545) in which the parameterization of near bed wave orbital velocities were modified according to Boechat . This modification improved the shallow wave-induced sediment transport and the resulting coastal morphology. Most importantly, this improvement allowed the modelling of coastal sediment transport with a more realistic coupling between the cross-shore and alongshore sediment transport fluxes and consequently a better overall sediment balance (see Boechat Albernaz et al., 2019). The Delft3D model is an extensively applied morphodynamic model for finite difference solving of the momentum and continuity equations for unsteady flow in depth-averaged mode through the Navier-Stokes equation with hydrostatic pressure approximation (Deltares, 2017). Waves were simulated with embedded Delft3D-WAVE (i.e., SWAN spectral wave model, see Booij et al., 1999;Ris et al., 1999) online coupled (i.e., two-way interaction) to the FLOW module which computes wave-current hydrodynamics, sediment transport, effects of vegetation, and morphological evolution (Deltares, 2017;Lesser et al., 2004).
We applied the VR04 sediment transport predictor ( van Rijn, 2007avan Rijn, , 2007bvan Rijn et al., 2004, improved according to Boechat  for the computation of near-bed wave orbital velocities. We choose the VR04 predictor because it is well-calibrated for wave-and current-driven sediment transport, includes wave-current interaction and intrawave sediment transport, and conceptually separates bed and suspended load for currents and waves. The mud fractions are treated as cohesive sediments and the deposition and erosional fluxes are computed according to Partheniades-Krone formulation (Partheniades, 1965) based on user-defined critical shear stresses for erosion and sedimentation. The critical shear stress for erosion was set to 0.5 N/m 2 and the sedimentation threshold to 1000 N/m 2 (the high value means that it always allows for sedimentation), both default values. For simulating the eco-morphodynamic effects of vegetation (eco-engineering species), we online coupled (i.e., two-way interaction) Delft3D to a dynamic vegetation module adapted from Brückner et al. (2019) that computes vegetation settlement and mortality derived from the hydro-morphodynamic development. In turn, the vegetation alters the hydro-morphodynamics through vegetation density-dependent roughness and drag (Baptist et al., 2007), which completes the feedback loop. The model was inspired by the sediment-dominated coastal systems of the North Sea, where the accumulation of organic material can be neglected. In view of the model complexity and the imposed spatial limit of the tidal basin, the succession of vegetation and peat builders (i.e., found further landward) was neglected.
Below, we detail the relevant model settings, including the initial and boundary conditions, and the dynamic vegetation characteristics.

Model Domain
The model domain was broadly inspired by the East Frisian Wadden Sea (see Figure 1 and Fitzgerald et al., 1984) and consists of a 20 by 5.5 km idealized coast connected to a basin of 1 m depth and 4.4 by 7.7 km dimensions via a single 1 km wide inlet enclosed by coastal barriers (Figure 2). These dimensions, specially the shallow basin with 1 m depth, represents a developed state of a backbarrier basin filled with clastic sediments, where intertidal areas can form and vegetation can settle. In order to perform simulation of a near-equilibrium condition between the coastal profile, tides, and wave climate, the coastal domain has a time-spatial averaged profile from the Dutch coast on the North Sea, and the tides and wave climate were schematized based on long-term measurements on the North Sea.
The coastal bathymetry was calculated from an alongshore uniform profile derived from a time-spatial averaged nearshore bed profile obtained from the Dutch JARKUS database (Rijkswaterstaat, 2017) on the Holland Coast near Leiden-Katwijk, the Netherlands. The profiles were measured between 1965 and 2010 along 31 alongshore distributed locations covering 15 km of coastline from the dunes up to approximately 15 m depth. The average profile has a slope of 1V:185H up to 2.3 m depth, steepening to 1V:63H up to the mean water line (as in Boechat Albernaz et al. (2019)). A comprehensive morpho-sedimentary description of the Dutch coast is presented in Short (1992), Wijnberg and Terwindt (1995), and van Rijn (1997).
We included four sediment fractions, that is, 200 and 125 μm sand, silt, and clay, in order to represent the different environments from the wave-exposed coast to the protected and vegetated mud flats ( Figure 1). Each sediment fraction was individually recorded by the stratigraphic bed module from Van Kessel et al. (2012) to allow different sediment mixtures and the effect of differential bed composition on sediment transport rates for each sediment fraction. The module tracks and stores the bed composition with a user-defined vertical resolution, here 0.1 m, and the sediment transport is computed for the active top layer on the basis of the top sediment mixture. The sediment supply was defined at the open sea boundaries as equilibrium concentration supply for the sand fractions, while clay was supplied as user-defined concentrations varying between 15 and 60 mg/L (within measured range of, e.g., van Kessel et al. (2011) and Kleinhans et al. (2005)) according to the modelling scenarios (see Figure 3). Note, however, that the open sea boundaries are far away from the tidal inlet, so that the sand concentrations at the inlet are entirely determined by the dynamic interactions within the littoral zone.

Tides, Wave Climate, and Sea Level Rise
The tidal boundary conditions were inspired by the North Sea records. The model has alongshore propagating morphological tides (Cayocca, 2001;Latteux, 1995;Lesser, 2009) with 0.75 m M2 amplitude, 0.075 m M4, 0.035 m M6 components and a generic diurnal (i.e., 24 hours) component of 0.2 m to create spring-neap cycles. All tidal components were rounded to integer hours.
The wave conditions of the Dutch Coast were recorded at, and summarized for, the IJmuiden "Munitiestortplaats" directional buoy in the North Sea between 1990 and 2016 (Boechat Albernaz et al., 2019). The wave directions were rotated 30° to realign the waves with our idealized coast based on the local shoreline orientation (Figure 2d). The wave climate consist of short-period waves averaging 4.6 with 11.7 s maximum. The maximum wave height  Figure 2e). The wave direction has two main components from SW (200°) and NW (310°) but also has a significant frequency of parallel and offshore going waves for 13% of the year.
The wave record was reduced and schematized into 13 representative wave conditions to optimize the computational effort following Walstra et al. (2013) and Benedet et al. (2016). The wave reduction consisted of 4 directions and 3 heights plus an average wave condition for the duration of offshore-directed and shore-parallel waves. The representative wave conditions were obtained by dividing the wave climate into 12 bins of equal energy E, being ∼ 2 . The wave period, direction, and the total duration were calculated as the average within each energy bin. Consequently, one wave condition energetically represents a range, consisting of significant wave height, wave period, wave direction, and duration. The reduced wave climate is shown in Table 1, with condition number 7 being the averaged condition to replace the recorded offshore and parallel going waves.
The wave conditions were implemented as 13 independent model cases which were run in parallel while sharing the bed level every time step (12 s). At each merge event, the new bathymetry was calculated as the sum of morphological changes across the 13 models. To weigh the differing occurrence of each wave condition, we multiplied the morphological change by individual morphological acceleration factors (i.e., morfac, see Ranasinghe et al., 2011;Roelvink, 2006) from Table 1 that were derived from the duration of each wave bin. This approach means that less energetic but persistent wave conditions have a higher morfac while the more energetic (stormy) conditions are performed with a lower morfac, which represents the sporadic stormy conditions. This parallel mode, known as "mormerge" (called "parallel online approach" in Roelvink (2006)), allows for the usage of more parallel computational resources while also improving the model stability (Roelvink, 2006). This parallel technique also mitigates artificial effects of ordering the wave conditions that can affect the dynamics of bars (Walstra et al., 2012(Walstra et al., , 2013, which is important for the sediment bypassing and for the formation of channels and ebb delta sandy bars (Elias et al., 2019;FitzGerald et al., 2000). To further enhance model stability and avoid bias towards one tidal condition, the tidal signal of each case was shifted in phase such that the ebb and flood conditions with higher morfacs happen simultaneously in different model cases, and therefore minimize the net change within the tidal cycle (Roelvink, 2006). With this model approach, we capture the net long-term morphological evolution, rather than episodic events, with more efficiency for the timescale of one century.
Sea level rise commenced after 20 years of morphological development to allow for a similar incipient tidal basin, ebb delta, and vegetation pattern for all scenarios. Sea level rise was implemented as a linear water level rise of 0.5 m per 100 years (i.e., 5 mm/yr), considering the total simulation time. This value was inspired by the Holocene evolution of tidal basins in the western Netherlands. These transgressive basins formed before 7500 BP under high sea level rise rates, migrated landward, and filled in from 5000 BP onward under a decelerating rate of relative sea level rise (RSLR) combined with abundant supply of sediment (Beets & van der Spek, 2000;de Haas et al., 2018;Meijles et al., 2018). The rate of RSLR during the phase of basin expansion was 5 mm/yr and higher and dropped to less than 2 mm/yr during the phase of basin infilling (Hijma & Cohen, 2019). The applied rate of 5 mm/yr thus marks the transition from basin expansion to basin infilling. In addition, the 5 mm/yr rate of RSLR is also in line with the IPCC's projected climate scenario under RCP4.5 (Representative Concentration Pathway), which predicts 0.52 m of sea level rise for the Dutch Wadden Sea area by the end of this century (Vermeersen et al., 2018). This value is bracketed by RCP2.6 and RCP8.4, which predict 0.43 and 0.84 m of global sea level rise by the year of 2100, respectively (Oppenheimer et al., 2019). However, higher values of RSLR in the order of 5 mm/yr or higher are expected locally due to the added effects of, for example, subsidence and extraction of water and gas from the subsoil (e.g., van Dobben et al., 2022;Minderhoud et al., 2020). Therefore our choice of the SLR rate represents the past transition from transgression in the Dutch coastal systems towards infilling (Vos, 2015) but is also in line with the future projections of relative sea level rise of, for example, the Netherlands which also accounts for subsidence.

Dynamic Vegetation Model
To investigate the effects of saltmarsh-type vegetation on the morphodynamics and related infilling processes, we implemented in the model a generic saltmarsh species, similar to Brückner et al. (2019), that represents common saltmarsh compositions in the northwestern Europe. At each ecological time step coupling (defined at every neap-spring tidal cycle, i.e., 24 hr), the results of the Delft3D hydro-morphodynamics (i.e., inundation period, flow velocity, and bed level changes) were fed into the dynamic vegetation model to calculate the new vegetation composition based on rules that determine the establishment and mortality. To correctly account for the inundation period and velocities, the vegetation module was coupled to the results of the case 7 (Table 1), which represents the average wave condition. The new vegetation coverage and properties were then applied to all cases. Below, we describe the main parameters, and more details on the dynamic vegetation model are provided in Brückner et al. (2019).
To translate the vegetation effects into the hydro-morphodynamic model, we applied the Baptist predictor (Baptist et al., 2007) to calculate the effects of vegetation on hydraulic roughness and flow drag force. This combination allows for a more realistic sediment transport computation when compared to roughness-only methods which leads to the overprediction of sediment transport rates due to the increase of the bed shear stress.
The dynamic vegetation model included settling and mortality rules conceptually similar to those in the literature (Brückner et al., 2019;D'Alpaos, 2011;Fagherazzi et al., 2012;Kirwan & Murray, 2007;Marani et al., 2013;Oorschot et al., 2016) for mature vegetation. The interaction between vegetation and waves was not accounted for in our models as the wave action inside the basin, where vegetation is present, is negligible. Due to the long timescale of the model, the vegetation growth phase was neglected as the ecological time step in morphological time is larger than the plant growth timescale. Therefore, settling was applied for mature plants with a specific plant height of 2 m as emergent vegetation, a root length of 1 m, a stem diameter of 2 cm, and a stem density of 400 stems/m 2 (Bouma et al., 2013;Leonard et al., 2002;Yamasaki & Tange, 1981). The initial coverage was set to a fraction of 0.25 (i.e., at each settling event the vegetation can settle at 25% of the grid cell) and can reach up to 1.0 representing the maximum coverage at one grid cell. At each coupling, the vegetation can recolonize suitable numerical cells as long as the maximum total fraction of 1 in the cells is not exceeded. The vegetation fraction implicitly accounts for the effect of plant biomass in the hydrodynamics via the Baptist predictor. The mortality induced by flow velocity, desiccation, and inundation stress linearly reduces the vegetation fraction in a cell. In addition, bed level changes that exceed the plant size immediately remove the entire fraction representing either complete burial for sedimentation larger than the plant height or removal for erosion deeper than the root length.
Vegetation settles in grid cells with an initial fraction representing patch density rather than individual plant density and therefore implicitly accounts for rhizome growth of, for example, Spartina anglica and plant dispersion. New fractions can added in suitable cells at every coupling time step, which represents lateral expansion and seedling establishment and leads to saltmarsh growth such as observed in nature (Brückner et al., 2019;Fagherazzi et al., 2012;Marani et al., 2007;Schwarz et al., 2022). The mortality is a portion of the initial fraction, which allows for constant die-off rates independent of the plant fraction present. This strategy leads to dense vegetation higher up the marsh while lower densities populate the more stressful locations (Brückner et al., 2019). The vegetation model does not include organic accumulation, which is, in the sediment-dominant systems considered here, only important landward of the tidal basins covered in the model. Here, we considered mortality values typical for S. anglica, Phragmites australis, and Scirpus maritimus (Bouma et al., 2013;Leonard et al., 2002;Yamasaki & Tange, 1981). Succession of species at higher elevations was ignored. The causes of mortality were implemented as dose-effect relationships (Brückner et al., 2019)

Tidal Basin Biomorphodynamics
Our results show strikingly different morphodynamic end-members between the vegetated and unvegetated scenarios (Figure 3). The reference runs without SLR (models 1 and 5 in Figure 3) demonstrate that vegetation 10.1029/2022JF006703 9 of 22 increases the number of channels and distributaries inside the basin. Also, without vegetation the northern coastal barrier migrates southward (alongshore) 1.3 km (i.e., 13 m/yr migration rate) and confines the inlet, while the scenario with vegetation develops a more stable and wider inlet (Figure 4). In the absence of vegetation, a single main channel develops inside the basin and migrates southward accompanying the inlet migration, whereas vegetation creates more equally distributed and stable channels landward of the inlet.
All scenarios initially evolve with reducing tidal prism through the inlet (Figure 5a). However, while unvegetated basins continue reducing the tidal prism while importing sediments, the vegetated scenarios revert to a trend of increasing tidal prism, enlargement of the inlet cross-section ( Figure 5b) and export of sediments ( Figure 5c). This sediment export trend causes the growth of the ebb-tidal delta (Figure 5d) while the vegetation cover declines for the sea level rise scenarios (Figure 5e) as the vegetation cannot keep up with the increasing environmental pressure. The vegetation coverage was only stable with constant sea level. The accommodation space, calculated from the volume of water inside the basin, stays constant for most scenarios (Figure 5f), except for the combinations of sea level rise and vegetation (models 6-8). The constant and stable accommodation space indicates a steady state condition between the creation of space driven by the sea level rise and the sediment infilling. The combination of sea level rise with various mud concentrations does not modify the general evolution of the unvegetated scenarios (models 2-4) when the basin imports more sediments to keep up with the increase of accommodation space. However, the corresponding scenarios with vegetation (models 6-8) show large deviations from the reference scenario (model 5). The accommodation space and tidal prism evolution of the vegetated reference scenario (model 5) is rather stable and of similar pattern as the unvegetated scenarios. The sea level rise in combination with vegetation results in underfilled basins and subsequent vegetation mortality in the late evolution stages. The increase in mud concentration, for example, in scenario 8 with 60 mg/L, results in more shoals and vegetation patches (see Figure 3) in comparison with scenario 6 and 7 (with 15 and 30 mg/L of mud, respectively). Nonetheless, the increase of mud supply is not sufficient to counteract the general vegetation collapse and drowning trend.
The vegetation in the basin triggers a cascade of positive feedbacks that lead to net erosion followed by vegetation mortality and basin drowning. The cascading feedbacks develop as follows: (a) an increase of the average depth  of the shoals and channels (Figure 6), as shoals did not keep up from the beginning of SLR (25-30 years), while the main channels got deeper after 40 years (due to vegetation and due to increasing tidal prism, see below), (b) a shift from sediment import to export around 50-60 years as the increase of flood storage led to ebb-dominance ( Figures 5 and 7), and (c) sediment starvation and erosion of the distal (landward) reaches of the basin, despite connection to tidal channels, which resulted in vegetation mortality and the formation of deepening ponding areas around 80 years (Figures 3 and 7). Conversely, the basin without vegetation continuously imports sediments and grows higher shoals.
The long-term erosional trend associated with sea level rise in vegetated basins (Figure 7) suggests that the local eco-engineering effects of saltmarsh plants can lead to system-scale drowning. Initially, sparse vegetation allows for local sediment trapping, with a typical pattern of onshore-directed transport on the tidal flats and further landward in the channels and offshore-directed transport in the deeper channels (Figure 7b at 20 years). The larger vegetation patches reduce sedimentation on the distal tidal flats, compared to the unvegetated cases. Moreover, the increased friction on shoals and intertidal areas combined with deeper channels increase flood storage and cause a change to ebb-dominance and consequently net sediment export and drowning under sea level rise. Important to note that the combination of vegetation and SLR results in the lack of steady state. Without SLR the vegetated basin reaches a steady state even with relatively low mud supply (i.e., 15 mg/L) and without vegetation the basin infills under SLR irrespective of the offshore mud supply. One could argue that scenario 5 (vegetation + no SLR + [15] mg/L) also shows a sediment export trend, similar to the vegetated scenarios under SLR that lead to drowning. However, the export of sediments in scenario 5 strongly decreases its rate after 80 years, which indicates a convergence to a steady state condition. Furthermore, the accommodation space and tidal prism of scenario 5 are rather stable and of similar trend and magnitude of the unvegetated cases that likewise reached steady state. On the other hand, the sharp increase in both tidal prism and accommodation space for the scenarios 6-8 indicate that the barrier tidal basin is not in steady state when combining SLR and vegetation.

Transition From Infilling to Drowning
The question is what controls the system-scale infilling or drowning of the tidal basins under the imposed sea level rise. The change from import (flood-dominance) to export (ebb-dominance) with ongoing sea level rise takes place before the vegetation starts to collapse. Figure 8a shows that first the tidal basin changes from importing to exporting sediment, induced by the eco-engineering effect of vegetation and later the vegetation coverage starts to decrease due to erosion (uprooting) and drowning (longer inundation periods). In other words, vegetation mortality results from initial erosion then later drowning. The only exception is with low mud supply (15 mg/L) when the erosion and vegetation mortality happen simultaneously.
The presence of vegetation reduces the basin dynamics (Figure 8b), defined as tidal prism (i.e., tidal discharge through the inlet) divided by accommodation space (i.e., total water volume in the basin), when compared to the nonvegetated analog scenarios. That means a decrease of tidal prism relative to the accommodation space which translates in more water retention, larger flood storage, and therefore ebb dominance. However, the nonvegetated scenarios show a persistent infilling trend, whereby the higher the mud content, the higher the infilling rate. Conversely, vegetation effects reduce the infilling independently of mud concentration. Figure 8b shows that scenario 8 (vegetation + SLR + [60] mg/L) has an incipient infilling trend but reverts to the same lower infilling rate as the other scenarios. In summary, from our set of models we observe a cascade of biomorphodynamic events that in the long-term resulted in distinct end-members (Figure 3) of the coastal barrier system following basin development and vegetation settling with sea level rise. The unvegetated basins were able to infill with sediments under all scenarios of sediment supply and sea level rise towards a steady state regarding tidal prism and accommodation space. Conversely, our vegetated basins only achieved a steady state without sea level rise. The increase of accommodation space induced by sea level rise followed an increase in tidal prism that could not be compensated by sedimentation. In fact, the basin developed ebb-dominance which resulted in net erosion and vegetation mortality, which represents a positive feedback towards erosion and drowning.

Discussion
Our novel set of models show that vegetation changes the morphodynamic feedbacks of the coastal barrier and that under sea level rise the basin drowns and the marsh coverage declines irrespective of the sediment supply.  Figure 3. Now, we discuss the different (quasi-) steady states, how our model results relate to past and modern natural environments, and how our results enlighten future adaptations in face of climate change and sea level rise.

Multiscale Biomorphodynamic Interactions and Steady States
Previous works showed that vegetation dissipates hydrodynamic energy (e.g., van Wesenbeeck et al., 2022) and consequently increases sedimentation (e.g., Contti Neto et al., 2022;D'Alpaos et al., 2007;Fagherazzi et al., 2012;Nardin et al., 2018). This mechanism is supposed to stabilise landscapes and to allow for vertical bed level accretion to keep up with sea level rise , suggesting that restoration of wetlands is an invariable nature-based solution for coastal protection (Schuerch et al., 2018). This contrasts with the main findings of this study. Our results show that the long-term biomorphodynamic interaction of vegetation and sea level rise may lead to net erosion, vegetation mortality, and final basin drowning, while without vegetation the basin may keep up with sea level rise.
Vegetation initially reduces the basin dynamics (similar to Boechat Albernaz et al. (2020)) and sedimentation takes place near feeder channels, which suggests that vegetation promotes fast infilling and reduction of accommodation space. However, similarly to D'Alpaos et al. (2007), Boechat Albernaz et al. (2020), and van Dobben et al. (2022), we showed that large patches of vegetation limit the conveyance of sediments, including mud, further from the feeder channels. This results in vegetation confining flow and sediment transport to channels, which was also observed in coastal (Schwarz et al., 2018;van Dobben et al., 2022) and fluvial systems . However, our models show that these local effects of vegetation have a basin-scale implication in the long-term. The long-term and large-scale morphodynamic effects of vegetation are to change the tidal dynamics through depth distribution ( Figure 6a) and flood storage within the basin (Friedrichs & Aubrey, 1988) and consequently the sediment balance (Figures 5c and 7b) of the entire basin.
Previous studies found that the tipping point for marsh survival and basin infilling are related to ratios of the sediment supply and sea level rise rates, for example, Kirwan and Murray (2007), Kirwan et al. (2010), D'Alpaos (2011), and Fagherazzi et al. (2012). However, unlike previous saltmarsh models, with enforced sediment concentration at the inlet boundary (e.g., Kirwan & Murray, 2007;Mariotti & Canestrelli, 2017), our model has a dynamic and free-evolving coastal zone where the wave-tidal conditions interacting with the nearshore and ebb-tidal delta to determine the sand-mud concentrations and the sediment budget (including sediment bypassing) of the basin. As a result, we were able to unravel a new and system-scale biomorphodynamic feedback where vegetation combined with SLR triggers a cascade of effects towards drowning that could not be counteracted by increasing sediment (mud) availability. These feedbacks eventually led to erosion and basin drowning under sea level rise even before the vegetation started to decline. Moreover, the basin drowning could not be counteracted by increasing the supply of mud as the basin develops a net sediment export behaviour. The long-term erosive trend resulted in the growth of subtidal flats, also called "pools" (see Figure 1 and Mariotti, 2020) as a consequence of the habitat loss. The unvegetated cases with the same conditions developed to a steady state between accommodation space and infilling. These findings demonstrate that the basin infilling and steady state is a result from the multiscale biomorphodynamic interaction related to the spatial connectivity of sediment pathways within the barrier coast (Elias et al., 2019;Pearson et al., 2020). In other words, the dynamics between the coastal area and the basin biomorphodynamics are determinant for the tidal basin and marsh survival.
This leads to the question whether these alternative equilibrium states are observed in natural systems. Figure 9 shows the depth distributions of natural tidal basins and modeled tidal basins, subdivided in vegetated and unvegetated conditions. From the natural systems, we observe the sharp predominance of inter-subtidal flats (leptokurtic distribution) and shallower wider channels on the unvegetated German and Dutch (Ameland) Wadden Sea (Figure 1 and Elias et al., 2019;Fitzgerald et al., 1984). Conversely, the vegetated systems (see Carrasco et al., 2018;Donatelli et al., 2020;Hein et al., 2012;Nardin et al., 2018), in general, present a wider depth distribution (platykurtic) or even bimodal distribution (e.g., see VCR in Figure 1) of deeper and narrower channels and a wider depth range between supra and subtidal flats including the higher (vegetated) shoals. This comparison strongly suggests that vegetation affects the depth distribution of channels, shoals, and flats, and that the alternative stable states exist in nature, similar to Marani et al. (2013). This does not exclude the possibility that other causes determine system development, equilibrium conditions and the presence and absence of vegeta tion. Nonetheless, our models support that the existence of different equilibrium states may manifest solely upon the presence-absence of vegetation while all other boundary and initial conditions are the same.
The tidal basin evolution and its response to changes in conditions, such as SLR or marsh decline, also relates to the larger scale coastal barrier behavior. Over long timescales (O-centuries to millennia), a barrier coast may adjust the accommodation space of the tidal basin via the translation of the barrier (i.e., regression or transgression) (FitzGerald et al., 2018;van der Spek & Beets, 1992;Vos, 2015). The FitzGerald et al. conceptual model (see Figure 22 in FitzGerald et al. (2018)) shows that marsh decline follows an increase in the tidal prism, accommodation space, and the growth of the ebb-flood-deltas. However, at a later stage, the barrier tends to migrate landward (i.e., barrier transgression) which reduces the accommodation space. Our models, despite the fact that cannot simulate barrier migration, show similar response regarding the marsh collapse and the tidal prism growth. However, in our models the marsh decline results from tidal prism growth and increasing erosion that started decades before the vegetation decline.
The revealed multiple equilibrium states and complex multiscale feedback mechanisms controlling vegetated and unvegetated tidal basins connected to the barrier coast points at the risk of using simplified models or empirical equilibrium relations to predict their morphological evolution. The empirical equilibrium relations (e.g., inlet cross section vs. tidal prism) or models with boundary conditions imposed inside the barrier coast (e.g., at the inlet) are insufficient to predict long-term effects of sea level rise on barrier coasts. As we showed, the complex barrier coast multiscale feedbacks and sediment pathways are crucial in predicting the ebb-tidal delta and the inlet geometry and the resulting sediment fluxes across the inlet, which are indispensable to understand the system-scale response of barrier coasts to SLR and sediment supply. A comparison of our model results and existing classical relationships between tidal prism and cross-sectional area of the inlet (Figure 10) shows that the reference runs without sea level rise have a much steeper trend than the empirical relationships, and in fact do not fully converge on any of them. Moreover, the final situation (top inset in Figure 10) after a century differs much between the runs, with vegetated runs evolving towards larger cross-sectional inlet areas than unvegetated runs, even without sea level rise. The models without sea level rise, for which the simple equilibrium relationships were derived, have different stable states due to the presence or absence of vegetation, in agreement with observations (Marani et al., 2013;Schuerch et al., 2018;van Belzen et al., 2017). To reveal those systemic connections across the subdomains of a barrier coast the modelling of the entire barrier coast is required.

Past and Future Evolution of Barrier Coasts
The model results of basin evolution under SLR and the biomorphodynamic feedbacks are consistent with the paleorecords of the mid-Holocene Dutch coastal plain that inspired this research. In our models with sea level rise on a late stage of infilling (i.e., mature infilled basin), the vegetation perished and lost its infilling effects for a sea level rise rate of 0.5 m per century, despite the littoral (sand) and offshore (mud) supply of sediment. The vegetated basin under SLR could not attain a steady state condition with respect to tidal prism and accommodation space (Figures 5 and 8b). Such (quasi-) steady state condition was only attained without sea level rise for both vegetated and unvegetated scenarios and for the unvegetated scenarios with sea level rise. We expect that lower SLR rates would show similar behaviour to our "no SLR" scenarios. Along the present-day Holland Coast, unvegetated backbarrier basins formed between the present-day cities of The Hague and Alkmaar around 7500 yr BP. These barrier systems formed under a sea level rise rate of approximately 0.3-0.5 m per century, while vegetation began to expand around 5000 yr BP when the sea level rise rate declined to 0.1-0.05 m per century Hijma & Cohen, 2019;Pierik et al., 2022;Vos, 2015). The paleogeographical reconstruction of the basin around Leiden-Katwijk in the Old Rhine area from Pierik et al. (2022) showed that vegetation settled on the margins of the tidal basin under higher sea level rise rates and expanded towards the mouth as the rate of sea level rise declined. The decline in RSLR rate together with the substantial fluvial sediment supply from the Rhine river promoted the infilling of the basin that allowed for the settlement and expansion of vegetation.
Regardless of predicting a precise SLR rate for the vegetation collapse (which is not our intention with the schematised model setup and the few runs with limited SLR scenarios), our results highlight the possible fate of coastal barrier systems under sea level rise, which is projected at 0.43 m rise until 2100 for RCP2.6 and 0.84 m until 2100 for RCP8.4 (Oppenheimer et al., 2019). Considering the predicted range of SLR, the reconstruction of the past Holocene Dutch coast and modern marsh research (FitzGerald et al., 2021) suggests that supratidal and intertidal areas will likely be reduced under RSLR rates in the order of 5 mm/yr. In this scenario, the vegetation coverage will also decline due to drowning, similar to the model results. The loss of supratidal and intertidal areas and vegetation is especially valid for embanked and diked areas where the retreat and the landward expansion of the natural system is limited or impossible. If the past SLR rates and sedimentation were indeed the key limiting factor for vegetation expansion, then a strategy of nature restoration for increasing land-surface elevations may only work in sediment-rich environments, such as river estuaries (see Pierik et al., 2022), while tidal basins may not be able to counteract drowning induced by RSLR. It is important to note that organic accretion promoted by vegetation, not incorporated in our models, but present in other model studies, for example, (Fagherazzi et al., 2012), peat and soil formation is further subject to higher compaction rates in comparison with sand and mud. The higher soil compaction rate of organic-rich layer is mainly induced by the weight of clastic sediments, organic decay, and also from drying out when the water table is lowered by natural or human processes (Minderhoud et al., 2020;Pierik et al., 2018). That said, vegetation and organic accretion can be seen as a comparatively fast short-term mechanism of infilling (Lorenzo-Trueba et al., 2012) but perhaps not a long-term solution for land rising. Conversely, sand and mud clastic deposits are less subjected to compaction and less prone to be altered by decomposition and decay in comparison with the organic matter.
The collapse of vegetation observed in our models, irrespective of varying sediment supply, and the indications of an initiation of collapse in field sites (e.g., Baptist et al., 2016;FitzGerald et al., 2021;van Dobben et al., 2022) raise the question of whether the modelled saltmarsh vulnerability to sea level rise is oversensitive to model assumptions and short-term field measurements or if it is genuinely a fragile ecosystem that would indeed perish under accelerated sea level rise conditions (Kirwan & Megonigal, 2013). For example, our models do not incorporate organic soil growth or multiple marsh species  that would likely increase the marsh resilience and the ability of sedimentation and bed level rise (FitzGerald et al., 2018). However, our results can be seen as conservative estimates to be compared to sediment-dominated natural systems. Yet, this type of collapse can hardly be explored from paleogeographical records as there are no available record of eroded systems, let alone one where vegetation collapse was followed by erosion. Despite the lack of long-term detailed paleorecords, the collapse of vegetation and its effects have been observed and projected in natural systems for the coming decades. For example, at the well-studied marsh of Plum Island in the USA, the marsh area is forecast to reduce dramatically by 2050 (FitzGerald et al., 2021). The projection of marsh reduction in FitzGerald et al. (2021) points that the Plum Island marsh (dominated by organic accretion) may lose a great portion of its high marsh in the coming decades under the RCP4.5 scenario. Part of the high marsh will possibly turn into low intertidal marsh which in turn increases the tidal prism that, similar to our models, could drive the whole marsh system towards a collapse. Our models show that the increase of flood storage and tidal prism can be an earlier trigger for erosion and marsh decline. Likewise, the Venice Lagoon in Italy shows a complex distribution of elevation-dependent saltmarsh states (Marani et al., 2013), all of which collapse under sea level rise rates between 0.39 and 0.59 m per century (Carniello et al., 2009;Marani et al., 2007Marani et al., , 2013. In contrast, unvegetated systems (or mostly unvegetated) such as the Dutch Wadden Sea are recognized to be able to cope with SLR rates up to 4-5 mm/yr due to a net import of sediments (Huismans et al., 2022;van Dobben et al., 2022).
The evolution of tidal basins is intricately connected to other elements of the barrier coast such as the adjacent sandy coasts, inlet, and ebb delta (Elias et al., 2019;FitzGerald et al., 2000FitzGerald et al., , 2021Lenstra et al., 2019;Pearson et al., 2020;Robbins et al., 2022;van der Spek & Beets, 1992). The underfilled state of our vegetated tidal basins in comparison with the vegetation-bare basins means that vegetated tidal basins need less sediments to reach (quasi-) steady state. Therefore, more sediments, especially sand, are available to the adjacent coast and to the ebb-tidal delta. The vegetated scenarios after 100 years promoted larger ebb-tidal deltas in the order of 10 million cubic meters of sand and fed the downdrift coast with 20 million cubic meters of sediments in comparison with the unvegetated scenarios. Similarly, the main channels in our vegetated basins are approximately 0.5 m deeper on average (Figures 6 and 5). These changes promoted by vegetation mean better navigability for harbors and more sustainable downdrift coasts in terms of sediment availability. The lower trapping of sediments in the basin, which no longer acts as a sink of sediments, promotes more sediment bypass to the downdrift coasts. These eco-engineering effects of marshes, therefore, translate into societal benefits within the basin but also beyond the barrier coast domain. These benefits enclose, for example, less need of extensive dredging, mitigation of downdrift structural coastal erosion that minimize the need of beach nourishments to try and avoid the coastal squeeze and having to abandon coastal plains (Siders et al., 2019). Therefore, apart from the intangible value of vegetation, the biomorphodynamic effects of a healthy tidal basin have positive implications on the barrier system as a whole. Yet, the collapse of coastal vegetation, such as saltmarshes and mangroves, due to accelerated sea level rise and climate change, may not only trigger losses in important ecosystem services, but also changes in the morphodynamic equilibrium state towards an unknown alternative state that mostly likely differ from the analog unvegetated state due to hysteresis .

Conclusions
Our set of biomorphodynamic models of an entire barrier coast system with a vegetated backbarrier tidal basin show that vegetation drives the tidal basin morphodynamics to an alternative steady state condition in comparison with unvegetated scenarios under otherwise the same conditions. The local vegetation effects on hydrodynamics and sediment transport first reduces shoal growth while deepening the channels, which later cascades to a basin-scale transition from flood-dominant to ebb-dominant conditions and consequently a tipping from net import to net export of sediments.
The contrasting morphodynamic evolution between unvegetated and vegetated basins has major implications in the scenarios with sea level rise. Here, the vegetated scenarios, despite abundant sediment (mud) supply, could not keep up with the increasing sea level. The export of sediments combined with increasing water levels leads to an increase in tidal prism, erosion, and vegetation mortality. Vegetation mortality and the formation of subtidal flats (ponds) causes a further increase of tidal prism as a positive feedback mechanism towards basin drowning. Large-scale vegetation mortality ensued and the system drastically changes to a new, drowned state. Conversely, the unvegetated scenarios persistently import sediments under sea level rise that compensates for the increase in accommodation space.
These findings highlight the potential consequences for natural systems that cannot keep up with the increasing climate and sea level rise pressure. While vegetation has been argued to invariably enhance sedimentation to keep up with sea level rise, our results suggest that in the long term, vegetation can reduce sedimentation in tidal basins. This counter-intuitive and novel insight derives from a cascade of biomorphodynamic effects and feedbacks that revert the initially intertidal importing system to an export intertidal and subtidal state.
From a societal-economical perspective, the vegetated systems result in deeper channels and need less sediments to achieve a steady state condition. This implies that more sediments are available to the coastal system. This means less need for dredging of navigation channels and also less need for nourishment of the downdrift coast. Moreover, marshes and coastal vegetation should be preserved but should not be seen as a direct measure for sedimentation that can be indiscriminately applied in future adaptations to sea level rise. Rather, human interventions will need to account for the complex biomorphodynamics that determine the fate of tidal basins and barrier coasts.

Data Availability Statement
Delft3D steering settings from our reference scenarios (model 1 and model 5) and main model results are available at the repository YODA (Boechat Albernaz, 2022). Delft3D source code is freely distributed and available at the Deltares (SVN) repository from Boechat Albernaz (2019). The vegetation module is also available at Brückner (2020) based on Brückner et al. (2019). Data from natural systems (see Figure 9) were obtained from DGT (2011), Richardson et al. (2018), Donatelli et al. (2020), and Sievers et al. (2020).