Aggregated morphodynamic modelling of tidal inlets and estuaries

Zheng Bing Wang,Ian Townend*,Marcel Stive

aFaculty of Civil Engineering and Geosciences,Delft University of Technology,Delft 2600 GA,the Netherlands

bCollege of Harbor,Coastal and Offshore Engineering,Hohai University,Nanjing 210098,China

cUnit of Marine and Coastal Systems,Deltares,Delft 2600 MH,the Netherlands

dSchool of Ocean and Earth Sciences,University of Southampton,Southampton SO14 3ZH,UK

Abstract

Keywords:Tidal inlet and estuary;Morphodynamic modelling;Aggregation;Temporal and spatial scales;Sea level rise;Sediment transport

1.Introduction

Morphodynamic modelling ofcoastaland estuarine systems has developed rapidly over the past 30 years,with many well-established modelling tools (e.g.,Delft3D,FVCOM,and ASMITA),accompanied by guidance on the application of such models(e.g.,Roelvink and Reniers,2012,www.estuary-guide.net).Within this discipline,most research has focused on the development of numerical models that solve the flow and sediment transport equations to predict the change in morphology over estuary space scales and time scales from years to decades(e.g.,Dam et al.,2016;van der Wegen et al.,2009).These models can be used to examine the evolution of tidal creeks,intertidal flats,and channels at the local scale,where the spatial resolution is on the order of tens of metres and the temporal resolution is seconds to hours.For morphological modelling,this is usually taken to be the micro-scale.Cowell et al.(2003)proposed a scale classification as follows:micro-scale(with a temporal scale from one second to one year and a spatial scale from 1 m to 1 km),meso-scale(with a temporal scale from one year to one thousand years and a spatial scale from 1 m up to 1 000 km),and macro-scale(with a temporal scale of millennia and a spatial scale greater than 1 000 km).

The models used at the micro-scale are typically governed by fast hydraulic flow simulations,which adopt some forms of the empirical relationship to deal with sub-scale processes,such as turbulence and bed friction.However,this is just one level of abstraction.Another well-established class of models is concerned with the meso-scale,where the focus is on the state of morphological features or the morphological system.This can be both qualitative for describing the state and system behaviour and quantitative for providing predictions at the feature/system scale.For the latter,the resolutions of space and time are usually on the order of kilometres and decades to centuries,respectively.

The so-called regime theory is an example of analysis at the meso-scale.This has been extensively applied in rivers,estuaries,and deltas to describe key dimensions(e.g.,width,depth,cross-sectional area,and ebb-tidal delta volume)based on the hydraulics and in particular the characteristic discharge of rivers or tidal channels(Lacey,1930;Spearman,2007).In the 1950s and 1960s,energy and entropy arguments were proposed as a theoretical basis to define the coefficients used in the earlier empirical hydraulic regime relationships,although this remains an area of some controversy(Langbein,1963;Leopold and Langbein,1962).Most recently,methods based on equilibrium sediment transport conditions have been used to define the system state(or regime)that is in equilibrium(Di Silvio,1989;Stive and Wang,2003).Although these models solve process equations that are similar to those used in micro-scale models,as described above,the integration is over larger(meso)spatial and temporal scales,and they have been variously referred to as aggregated,box,and semi-empirical models in the literature.However,the classification of those models as semiempirical is misleading because all models across all scales use some forms of parametric closure based on empirical relations to solve the governing equations.The closure of the governing equations is generally at an appropriate level of the process formulation,e.g.,for small-scale models,the closure is at the level of the friction and/or the turbulence;for meso-scale models,the closure is at the level of the aggregated sediment transport;and for large-scale models,the closure is at the level of the morphological elements,e.g.,the entire ebb-tidal delta of a tidal inlet.

For morphological modelling,the meso-scale class of models has proved particularly useful when examining a longterm change related to engineering and coastal management problems.However,there remains a lack of clarity about how this type of models relates to models used at the micro-scale.To address this,we provide a detailed account of the progressive stages of aggregation used to move from a three-dimensional(3D)spatial description to an aggregated box model description.We then use the resulting mathematical description to illustrate the insights that aggregation at the macro-scale can offer.Whilst these insights could probably have also been developed at the micro-scale,this would require huge numbers of simulations and would still lack the clarity offered by the relatively simple expressions developed through an aggregation at an appropriate scale,i.e.,the scale of interest.

2.Stages of aggregation(from micro-to meso-scales)

In order to clarify the link between micro-scale and mesoscale models,we made use of the convection-diffusion equation,as used in many coastal and estuarine morphological models,to examine the implications of aggregation progressing from a 3D spatial description(which itself involves aggregation in the closure terms)to a model representing a system as a single box that exchanges sediment with the outside environment.The aggregation steps are described below.

2.1.3D to 2D aggregation at micro-scale

We start with the 3D convection-diffusion equation governing the sediment concentration:

where c is the sediment concentration;t is time;ux,uy,and uzare the flow velocity components in the x-,y-,and z-directions,respectively; εx, εy,and εzare the turbulent diffusion coefficients in the x-,y-,and z-directions,respectively;and wsis the settling velocity.This is in fact the mass conservation equation for sediment transport,forming the kinetic part of the theory for suspended sediment transport.The dynamic part of the theory for suspended sediment transport is reflected by the bed boundary condition.At the boundary near the bed,the sedimentconcentration,orthe sedimentconcentration gradient,can be prescribed.Their values can be calculated using a sediment transport formula(in the case of sand)or a formulation for the erosion rate(in the case of mud).

Eq.(1)can be integrated in the vertical direction to obtain the depth-averaged advection-diffusion equation:

where the over-bar represents depth-averaging of the corresponding variable;h is the water depth; αxand αyare the coefficients accounting for the effects of the flow velocity and sediment concentration distributions in the x-and y-directions,respectively;Dxand Dyare the dispersion coefficients in the xand y-directions,respectively;and fbis the sediment exchange flux between the bed and the water column.

This can be considered the first level of aggregation used by many micro-scale models.Note that the diffusion coefficients in the horizontal direction become dispersion coefficients as they also represent the mixing due to the non-uniform vertical distributions of flow velocity and sediment concentration.

2.2.2D to 1D aggregation at micro-scale

The next level of aggregation can be carried out by integrating Eq.(2)over the width of a river,estuary,or channel:

where the over-bar represents averaging over the cross-section,A is the cross-sectional area,α is the coefficient accounting for the effects of the distributions of flow velocity and sediment concentration within the cross-section,and Fbis obtained through integration of fbover the width.This is the equation for suspended sediment concentration used in 1D network models.

Eq.(3)can also be written as

where sxand syare the suspended sediment transport rates in the x-and y-directions,respectively.This equation can be further aggregated by integrating over a part of or the whole estuary(a morphological element).Using Green's theorem,the integration yields

where V is the volume of the water body of the area,Jiis the sediment transport at the open boundary i(a positive sign means that the transport is directed to an element,i.e.,import),FBis obtained through integration of fbover a morphological element,and C is the average sediment concentration in the element.This equation can also be derived by considering the mass balance of sediment in the whole water body.

2.3.Furtheraggregationinspaceandtimeatmeso-scale

Eqs.(2)-(5)can also be aggregated in time,e.g.,over a tidal period or a much longer time.As an example,integration of Eq.(3)over a tidal period yields

The first term in Eq.(2),representing the change of sediment storage in the water column,is neglected as it becomes much less important at longer time scales than the term on the right-hand side representing the sediment exchange with the bed.All the other terms remain basically the same,but the parameters and variables now represent the tidally averaged values.The residual flow velocities(uxrand uyr)cause advection,and the tidal flow now becomes the major mixing agent for the dispersion represented by the coefficients Dxand Dy,as elaborated by Wang et al.(2008).

Aggregation of Eq.(5)in time yields the following equation used in an aggregated model for the aggregated scale morphological interaction between tidal basin and adjacent coast(ASMITA,Stive and Wang,2003):

2.4.Sediment exchange between bed and water column

The advection-diffusion equation describes a mass balance for suspended sediment transport,no matter what level of aggregation is adopted.It does not describe any dynamics but only the kinetics of suspended sediment transport,because the equation itself does not provide the information about the flux of sediment exchange between the bed and the water column.In its aggregated forms(Eqs.(2)-(7)),this flux is expressed with a term in the equation.It cannot be determined by solving the equation,but needs to be prescribed.In fact,this flux determines the morphological change as well,and its formulation represents the dynamic part of the suspended sediment transport model.

For the 3D form(Eq.(1)),the dynamic part is introduced via the bed boundary condition.At the boundary near the bed,the sediment concentration,the vertical gradient of the sediment concentration,or a linear combination of the two factors can be prescribed.When prescribing the sediment concentration,it is often argued that the concentration near the bed can instantaneously adjust to the local flow condition,and therefore the equilibrium value of sediment concentration can be prescribed.In fact,for steady uniform flow,the solution for the vertical profile of sediment concentration,far away from the boundary,is the equilibrium concentration profile.The prescribed value of sediment concentration at the bed boundary is therefore by definition the equilibrium concentration.A similar argument can be made for the prescribed concentration gradient at the bed.Physically,by prescribing the sediment concentration,the vertical sediment flux due to settling is given but the vertical flux due to turbulent mixing is left undefined.However,by prescribing the concentration gradient,the opposite is true.For non-cohesive sediment,the required equilibrium concentration or the equilibrium concentration gradient depends on the flow conditions and the sediment properties and can be determined with one of the many sediment transport formulas(e.g.,Soulsby,1997;van Rijn,1993).Note that the equilibrium values of sediment concentration and its gradient are related to each other as the downwards flux due to settling and the upwards flux due to turbulent mixing balance each other out at equilibrium.For cohesive sediment,deposition and erosion rates according to the Krone-Partheniades formulation are often used for the bed boundary condition.According to this formulation,deposition occurs only when bed shear stress τ is below a critical value τd,and erosion occurs only when bed shear stress is above another critical value τe.With this formulation it is always impossible to define the instantaneous equilibrium concentration by balancing the deposition and erosion rates.For τd＜τe,this leads to an equilibrium concentration equal to zero for τ＜ τd,an equilibrium concentration undetermined for τd＜τ＜τe,and an equilibrium concentration infinitely large for τ ＞ τe.However,if the hydrodynamic condition is fluctuating due to some factors,e.g.,tides,a tide-averaged equilibrium concentration can be defined.Recently,Winterwerp et al.(2012)argued that deposition always takes place,i.e.,τdshould be infinitely large.With this approach it is then possible to define the instantaneous equilibrium concentration for mud.

The advection-diffusion equations when aggregated in space but not in time(Eqs.(2)and(3))are also used in microscale models.The flux fbor Fbcan be calculated with the Krone-Partheniades formulation for cohesive sediment,or with a formulation derived from an asymptotic solution of the 3D advection-diffusion equation (Galappatti and Vreugdenhil,1985;Wang,1992)for non-cohesive sediment.These formulations require information concerning hydrodynamic conditions expressed with the flow velocity,bed shear stress τ,etc.,in addition to the sediment properties,e.g.,the median grain size of sediment D50:

When the advection-diffusion equations are aggregated in time,the detailed information on hydrodynamic conditions required by a sediment transport formula or the Krone-Partheniades formulation is no longer available.The sediment exchange flux between the bed and the water column needs to be expressed with available aggregated morphological and hydrodynamic parameters.This formulation,representing the dynamics of sediment,is dependent on the level of aggregation,which is a key difference from other models.The formulations often have the following form:

where γ is a dimensionless coefficient,and ceis the equilibrium sediment concentration.This is similar to the formulation of Galappatti and Vreugdenhil(1985)for micro-scale models.The difference between the models is the way in which ceis calculated.

Di Silvio(1989)related the equilibrium concentration to the morphological state variable.As an example,for shoals influenced by wind waves,he used

where n is a power.It is noted that this relation implicitly defines the morphological equilibrium,at which ceis equal to the sediment concentration prescribed at the open boundary of the model.

In the ESTMORF model(Wang et al.,1998)and the ASMITA model(Stive et al.,1998),the equilibrium concentration is dependent on the ratio between the morphological state variables and its equilibrium value.As an example,for the channels in the ESTMORF model,the formulation is

where Aeis the equilibrium value of the cross-sectional area,and ceis the global equilibrium concentration.When the whole system is in morphological equilibrium,c=ce=cE.The rationale behind this formulation is that the ratio Ae/A is an indication of the flow strength in the channel with respect to equilibrium,as Aeincreases with the increasing tidal prism.The formulation is thus analogous to a power law for sediment transport capacity.

3.Modelling at meso-scale

Describing the system in terms of larger-scale morphological elements has the advantage that empirical relations are available for morphological equilibrium.In contrast,the detailed micro-scale models often encounter difficulties in correctly reproducing,or even achieving,morphological equilibrium.However,the aggregation in the meso-scale schematisation also introduces difficulties in the description/representation of the hydrodynamics in sediment transport processes.To avoid such difficulties,the classic aggregated morphodynamic models assume exponential decays of disturbances with respect to morphological equilibrium(Eysink,1990).The decay time scales are often empirically derived.Therefore,these models have been called empirical models.To our knowledge,Di Silvio(1989)and Di Silvio et al.(2010)were the first to introduce a suspended sediment transport module into an aggregated(in time)morphodynamic model.The long-term averaged sediment concentration field is simulated directly with an advection-diffusion equation,just like a micro-scale model that solves the instantaneous sediment concentration field.In addition,the local equilibrium value for the long-term averaged suspended sediment concentration is defined and related to the local morphological state variable following some physical arguments,and again this is conceptually the same as a micro-scale model.As an example,the equilibrium value of the long-term averaged sediment concentration on shoals is argued to be inversely proportional to a power of the local water depth for given wind and wave climates.By relating the localequilibrium sedimentconcentration to the morphological equilibrium,the aggregated morphodynamic models for ESTMORF(Wang et al.,1998)and ASMITA(Spearman,2007;Stive et al.,1998;Stive and Wang,2003;Townend et al.,2016a)have been developed.These models haveproved to beeffectiveformodelling long-term morphological development in estuaries(Townend et al.,2007;Wang and Townend,2012)and tidal inlet systems(Kragtwijk et al.,2004;van Goor et al.,2003).They combine the empirical relations for morphological equilibrium with process modelling of suspended sediment transport and therefore have often been classified as hybrid or semiempirical models.

The previous section makes it clear that the essential difference between meso-scale and micro-scale models(e.g.,ASMITA and Delft3D)is the level of aggregation rather than the degree in which their formulations contain empirical elements.They both use the advection-diffusion equation for modelling the suspended sediment transport.The formulations concerning sediment dynamics,i.e.,the exchange between the bed and the water column,contain empirical elements in both models.The ASMITA model makes use of the empirical relations for morphological equilibrium and Delft3D uses a sediment transport formula,which is also empirical in character.The empirical elements in both models also involve the same degree of uncertainty.In this sense,it is not meaningful to classify the two types of models as respectively(semi-)empirical and process-based.These models only differ in the level of abstraction that is used to represent the real world.

In addition to the aggregation level and the formulation for the sediment exchange flux between the bed and the water column,the hydrodynamic module helps to differentiate various models.Obviously,the implemented hydrodynamic module depends to a large extent on the aggregation level of the model.As an example,the models for ESTMORF(Wang et al.,1998)and ASMITA(Stive et al.,1998)are based on the same type of formulations.However,the ESTMORF model can be coupled to a(micro-scale)1D network hydrodynamic model to simulate the required aggregated hydrodynamic parameters,e.g.,the tidal volume and tidal range,as various other hybrid models(e.g.,O'Connor et al.,1990;Wright and Townend,2006).In contrast,due to a higher level of spatial aggregation,the ASMITA model calculates the tidal prism for the prescribed tidal range and plan area.Micro-scale models can be 1D,2D,or 3D,depending on the hydrodynamic module.

The ASMITA formation was first proposed by Stive et al.(1998)for modelling long-term morphological development of tidal inlet systems in the Dutch Wadden Sea.A tidal inlet system is schematised into a limited number of morphological elements,at a level like an ebb-tidal delta,as described by Walton and Adams(1976).For each element,a water volume below a certain reference level or a sediment volume above a certain reference plane(not necessarily horizontal)acts as the integral state variable.A tidal inlet system is typically schematised into the following three elements(Fig.1):an ebb-tidal delta,with the total excess sediment volume relative to an undisturbed coastal bed profile(Vd)determined as its state variable;an intertidal flat area in the tidal basin,with the total sediment volume between the mean low water(MLW)and mean high water(MHW)(Vf)determined as its state variable;and a channel area in the tidal basin,with the total water volume below MLW(Vc)determined as its state variable.The hydrodynamic and morphological parameters are defined according to Lodder et al.(2019).

Fig.1 shows that the adjacent coastal areas,which can exchange sediment with the inlet system,are considered an open boundary.Tidal flats are depicted with a sediment volume(Fig.1(c)),but it is also possible to represent tidal flats with the volume of water over the tidal flats(Townend et al.,2016a).

Although commonly used for a tidal inlet system,this three-element schematisation is not the only possible schematisation.In fact,an ASMITA model can contain any number of inter-connected morphological elements.It is relevant to mention that the ASMITA model was applied to the Dutch Wadden Sea basins based on the assumption that the floodtidaldeltaextendsoverthe wholeback-barrierbasin,although it is not generally the case.The US east coast has many tidal basins,with a flood-tidal delta confined to a location near the inlet opening.This would require a different set of elements to represent the system.

Fig.1.Three-element schematisation for a tidal inlet system in ASMITA model.

Anyhow,a necessary requirement is that the morphological equilibrium of each element is defined and can be evaluated with the available(aggregated)hydrodynamic parameters in the model.In the following section,this three-element schematisation is used to explain the model formulation and to demonstrate some applications.In addition,a one-element model is used to explain certain points,in which a whole back-barrier basin of the tidal inlet system is considered a single morphological element with the water volume below MHW in the basin as the state variable.

The formulation of the ASMITA model has been detailed in various publications(Stive et al.,1998;van Goor et al.,2003;Kragtwijk et al.,2004;Townend et al.,2016a,2016b).For completeness,the basis of a three-element model suitable for the exploration of the impacts of sea level rise on tidal inlets and estuaries is provided in Appendix.

4.Impact of sea level rise on tidal inlets

Aggregation to a particular scale allows the modeller to focus on system states,or flows in and out of the system.Choosing an appropriate scale of aggregation can make certain problems more tractable,or provide greater insight into the system dependencies.This is illustrated by considering the sensitivity of tidal inlets to sea level rise,using the formulation derived for the meso-scale.For simplicity's sake,we adopted the ASMITA formulation for a single element(combined channel and tidal flat with no delta).Previous work examined the potential for tidal inlet systems to drown(Rossington et al.,2007;van Goor et al.,2003).The analysis that follows is similar but has more detail,to highlight how the various system parameters influence the sensitivity of the system to sea level rise.

The one-element ASMITA model for a tidal basin with the water volume below MHW,V,as the state variable,is as follows:

where Veis the equilibrium volume(without sea level rise)of the tidal basin,R is the sea level rise rate,Abis the horizontal area of the basin at MHW,w is the vertical sediment exchange velocity,and δ is the horizontal sediment exchange coefficient.The details about exchange and transport of sediment are explained in Appendix.A dynamic equilibrium is achieved for a constant rate of sea level rise if the volume V is constant,with dV/dt=0.

The dynamic equilibrium is possible only if R is smaller than a critical value Rc,with

For the dynamic equilibrium state,we have

This solution is depicted in Fig.2 for n=2.A linear model is derived by linearizing Eq.(12)around V=Ve,which is also shown in Fig.2:

where T is the morphological time scale,and

Fig.2.V/Vevs.R/Rc.

Comparison of Eqs.(13)and(16)yields a relation between T and Rcas follows:

A dynamic equilibrium for a constant rate of sea level rise is again achieved if V no longer changes with time,with dV/dt=0.This leads to

4.1.Assessingsensitivityofsystem torateofsealevelrise

Considering the system as a single element,a non-linear equation can be used to determine the sensitivity of the system to a given(constant)rate of sea level rise,as indicated by Eq.(13).The limiting rate of sea level rise is determined by four parameters,and it is instructive to examine how the limiting rate of sea level rise varies with these parameters.Typical ranges for each parameter are defined in Table 1.The way in which each of these parameters influences the limiting rate of sea level rise is shown in Fig.3.

The joint influence of the parameters is illustrated in Fig.4.Fig.4(a)shows the joint influence of the size of the basin and the availability of sediment.Fig.4(b)shows the joint influence of the horizontal sediment exchange coefficient and vertical sediment exchange velocity.

4.1.1.Relative importance of parameters

For the range of parameters considered,there is a similar degree of sensitivity of the basin size,the vertical sediment exchange velocity,and the horizontal sediment exchange coefficient to the limiting rate of sea level rise(Fig.3(a)-(c)).There is a much larger variation of the limiting rate of sea level rise for the range of equilibrium sediment concentration considered(Fig.3(d)).This range represents the state from a limited sediment supply to a relatively turbid system and highlights the importance of sediment supply.The basin size can have a significant effect on the system's ability to keep pace with the sea level rise.Large basins with a limited sediment supply are more likely to drown and small basins with a plentiful sediment supply should be able to maintain their form even if the rate of sea level rise increasesdramatically.In a similar manner,an increasing sediment supply increases the resilience of the system to a rapid sea level rise.Moreover,a large horizontal exchange coefficient and a high velocity of vertical sediment exchange with the bed may both reduce the vulnerability of the system to a rapid sea level rise.

Table 1Typical ranges for parameters that determine limiting rate of sea level rise for single element.

Fig.3.Plots to illustrate limiting rate of sea level rise varying with different parameters.

4.1.2.Sensitivity of system to sea level rise

Fjords,fjards,and rias all have limited sediment supply.Arguably,fjords and fjards are drowned systems.However,the potential consequence for rias,as found in the southwest of the UK,is that they can also be progressively drowned as a result of high rates of sea level rise.For coastal plain,tidal inlets,and drowned valley systems,the system setting is important.Systems under macro-tidal conditions are likely to be less vulnerable than those under micro-tidal conditions because of the higher horizontal sediment exchange rate.Similarly,sandy systems are likely to be less sensitive than muddy ones because the vertical sediment exchange rate is higher(sand has a larger settling velocity).However,in all cases,the availability of sediment is critical.The background supply may be suspended load or bed load and may be sourced from the river or the marine environment.In this model this is represented by the equilibrium sediment concentration,an important parameter for assessing whether the system remains healthy or degrades,as illustrated by the changes in Venice over the last two centuries(Townend,2010).

Fig.4.Limiting rate of sea level rise under joint influence of parameters.

It is important to recognize that Eq.(12)treats the estuarine system as a single element.This is,of course,a gross simplification.The ASMITA model subdivides the system into the main geomorphological elements,namely channel,tidalflat,and delta.Introducing this additional subdivision results in different responses for the different types of elements.In most cases,the tidal flats(comprising a smaller volume and being furthest from the source of sediment)are significantly more sensitive to the sea level rise than the tidal channel and delta components of the estuary.Consequently,several studies have shown that estuarine systems including tidal flats have the capacity to adapt to rates of sea level rise in the range of 5-30 mm/year(Rossington et al.,2007).

4.2.Application of ASMITA model

The ASMITA model is meant for simulating long-term large-scale morphological development of tidal inlet systems.It was originally developed for the tidal inlet systems in the Dutch Wadden Sea(Stive et al.,1998).The first applications to the Dutch Wadden Sea concerning the responses of the tidal inlet systems to human interferences(Kragtwijk et al.,2004)and to sea level rise(van Goor et al.,2003)were summarized by Stive and Wang(2003).Later,the model was applied to the Dutch Wadden Sea in environmental impact assessment studies for gas and salt mining and to other systems in the world,including estuaries(Townend et al.,2007;Rossington et al.,2007).It has also been implemented as a module in the DIVA model for evaluating the impact of climate change on the worldwide coastal system(Hinkel et al.,2013).More discussion on the ASMITA model and its applications can be found in Townend et al.(2016a,2016b).Here we discuss an application of the model to evaluating the impact of relative sea level rise on the Dutch Wadden Sea tidal inlet systems.

Fig.5.Tidal basins in Dutch Wadden Sea.

The Dutch Wadden Sea consists of a series of tidal inlets with their tidal basins(Fig.5).Multiple large-and small-scale human interventions in the past have basically fixed the boundaries of the Dutch Wadden Sea.The two recent major interventions were the closures of the Zuiderzee in 1932 and the Lauwerszee in 1969,and they are still influencing the morphological development(Elias et al.,2012;Wang et al.,2012,2018).Due to the closure of the Zuiderzee,the western part of the Dutch Wadden Sea(Texel,Eierlandse Gat,and Vlie inlets)is still far from morphodynamic equilibrium(Wang et al.,2018).

Historically,the morphological development of this system,consisting of barrier islands separated by tidal inlets and backbarrier basins with channels,tidal flats,and saltmarshes,has been closely related to the sea level rise.It was formed about 7 000 years ago under the influence of a rising sea level.The development of basins depends on the balance between sediment accretion and the relative sea level rise rate(Nichols,1989).Beyond a criticalrate of relative sea level rise at approximately 1 m per century,sediment import becomes insufficient,resulting in drowning of the basins(Carrasco et al.,2016;Deltacommissie,2008;Stive et al.,1990;van der Spek and Beets,1992;van Goor et al.,2003).The relative sea level rise is the sea level rise plus the sea floor subsidence due to,e.g.,gas and salt mining.

Anticipating that worldwide sea level rise will accelerate(KNMI,2017;Le Bars et al.,2017;Meehl et al.,2007),studies on the Dutch Wadden Sea have been carried out,resulting in a series of position papers coordinated by the Wadden Academy in the Netherlands.Vermeersen et al.(2018)projected three scenarios for future sea level rise in the Dutch Wadden Sea until 2100,following the representative concentration pathways(RCP2.6,RCP4.5,and RCP8.5).Fokker et al.(2018)predicted sea floor subsidence due to gas and salt mining in the coming decades.The morphological response of the Dutch Wadden Sea to the sea level rise and sea floor subsidence(together the relative sea level rise)was projected by Wang et al.(2018),on the basis of the ASMITA model results and supported by observations concerning the morphological development since 1926.The ASMITA model used was the same as the approach of van Goor et al.(2003),with the most up-to-date parameter settings of the ASMITA model for the tidal inlet systems of the Dutch Wadden Sea shown in Tables 2-4.

Following the approach of van Goor et al.(2003),the results are presented as the probability of drowning as a function of the rate of sea level rise,which has been obtained by considering all the relevant model parameters as stochastic variables.Based on this probability distribution,an uncertainty range for the critical rate of sea level rise can be obtained,as shown in Table 5.In this table,the value of the rate of sea level rise corresponding to a probability of 16% for drowning is the median value minus the standard deviation(if the distribution is a normal one).However,the values for the Texel and Vlie inlets appear to be unrealistic,as the lower limits are lower than the observed sedimentation rates in the past.Theobserved sedimentation rates at these two inlets are about the same as the sea level rise rate corresponding to the probability of 30% for drowning.Therefore,the range based on the probability of 30% is also given in Table 5.The median values of the critical sea level rise rate,corresponding to a probabilityof 50%,are used for the projections of the morphological response.

Table 2Tidal ranges and morphological parameters of three elements.

Table 3Parameters influencing morphological time scale.

Table 4Parameters for defining morphological equilibrium.

Table 5Critical sea level rise rates for tidal basins in Dutch Wadden Sea for five values of probabilities of drowning.

By comparing the median values in this table to the sea level rise scenarios presented by Vermeersen et al.(2018)and the sea floor subsidence projections by Fokker et al.(2018),Wang et al.(2018)concluded that the critical rate of relative sea level rise will not be exceeded for the RCP2.6 scenario in any of the basins in the Dutch Wadden Sea.For the RCP4.5 scenario,it will be exceeded only in the Vlie Inlet Basin starting from 2030;for the RCP8.5 scenario,it will be exceeded in the Texel Inlet Basin starting from 2050,in the Vlie Inlet Basin from 2030,and in the Amelander Inlet Basin from 2100.It is noted that exceedance of the critical value does not mean that a tidal basin will drown(loss of all intertidal flats)immediately.In fact,it is concluded that none of the considered tidal basins will drown before 2100,even for the RCP8.5 scenario.However,the accelerated sea level rise will cause loss of part of the intertidal flats.The loss will be noticeable in 2050 for the RCP8.5 scenario and in 2100 for the RCP4.5 scenario.

5.Conclusions

Aggregation is a valuable means by which a modeller can focus on representing the real world at an appropriate level of abstraction,choosing the attributes and system states of interest.A very crude abstraction may provide simplicity but insufficient resolution of the problem,whereas inclusion of too much detail may make the model computationally intensive,or even intractable.Hence the correct choice is important for efficiency but more important for the clarity of insight that can be obtained.

Using morphological models of tidal inlets and estuaries as an example,we have shown that the aggregation is a process of progressively grouping terms,or integrating over spatial and temporal scales.Thus,in the example presented,the models at the different scales(micro and meso)solve the same equations just at different levels of resolution.Importantly,some forms of empirical closure scheme are needed at all scales.Hence,the models are all empirical to some degree.

The insights to be obtained from aggregation to the scale of interest have been illustrated for the case of sea level rise.With some manipulation of the governing equations at the meso-scale it is possible to show how a tidal inlet subject to sea level rise is likely to respond,as a function of system size,tidal range,sediment supply,and type of sediment.Such a relatively simple examination of the key dimensions of the problem would be much harder to achieve using an equivalent micro-scale model.

Acknowledgements

This study was carried out within the framework of the Coastal Genesis 2 Project funded by Rijkswaterstaat,and cofunded by Netherlands Organization for Scientific Research(NWO)via the SEAWAD Project and the Dutch Royal Academy of Sciences via the project of Coping with Deltas in Transition within the framework of Programme Scientific Strategic Alliances between China and the Netherlands.

Appendix:Aggregated morphodynamic modelling of tidal inlets and estuaries

The bases of the hydrodynamics,morphological equilibrium state,sediment transport,and morphological changes used in the ASMITA model are detailed below.This tidal inlet sytem is presented in terms of three geomorphological elements,namely delta,channel,and tidal flat.The details of the formulation of an ASMITA model using multiple elements can be found in Kragtwijk et al.(2004)and Townend et al.(2016a,2016b).

A.1.Hydrodynamics

If the morphology of a tidal inlet system is only defined by the volumes of the three morphological elements mentioned above(see Fig.1),it does not need to model the hydrodynamics in detail.In fact,only the tidal range in the back-barrier basin and the tidal prism are needed to define the morphological equilibrium state.Tidal basins,such as those in the Dutch Wadden Sea,are relatively small with respect to the tidal wavelength.For such small basins it is acceptable to use the water volume between the mean low water(MLW)and mean high water(MHW)for the tidal prism.When the intertidal volume is defined in terms of sediment volume,the relation between the tidal range H and the tidal prism P is as follows:

where Abis the horizontal area of the basin at MHW,and Vfis the total sediment volume between MLWand MHW.The tidal range in the basin is dependent on the tidal range in the open sea and the morphological state of the tidal inlet system.The tidal range in the open sea should be considered a boundary condition.It is difficult to model the influence of the morphological changes on the tidal range due to the highly aggregated schematisation.However,for small basins like those in the Dutch Wadden Sea,the influence of the morphological changes on the tidal range is limited.Therefore,the tidal range in the basin H is prescribed.It is noted that it is still possible to have H varying with time,to take account of any trend or a cyclic change,e.g.,the nodal tide variation(Jeuken et al.,2003;Wang and Townend,2012).It is also noted that the model does consider the feedback from the morphological development to the hydrodynamics,via Eq.(A.1),because Vfis a morphological state variable.

A.2.Morphological equilibrium state

For all three elements,empirical relationships are required to define their morphological equilibrium dimensions.With these relationships the morphological equilibrium state of a tidal basin is fully determined if the basin area Aband the tidal range H are given.The equilibrium relationships can be prescribed in various ways(see Townend et al.,2016a,2016b).Here we outline the approach that has been used to study the dynamics of the Dutch Wadden Sea.

For the intertidal flat there are two empirical relationships used,one for its area and the other for its height(Renger and Partenscky,1974;Eysink and Biegel,1992):

where Afeis the equilibrium tidal flat surface area in m2,and hfeis the equilibrium height of the tidal flat.According to Eysink(1990),αfe= αf-0.24× 10-9Ab,with αf=0.41.

The equilibrium volume of the intertidal flat,i.e.,the equilibrium sediment volume between MLW and MHW,is given as

The channel volume Vcis defined as the water volume under MLW in the tidal basin.Its equilibrium value is related to the tidal prism as follows:

where αcis an empirical coefficient.

The equilibrium volume of the ebb-tidal delta is also related to the tidal prism(Walton and Adams,1976)as follows:

where αdis an empirical coefficient.

Using these equations with various empirical coefficients αcand αd,the morphological equilibrium of a tidal basin can be determined with two parameters,the total basin area Aband the tidal range H,with Ve=F（Ab,H）,where Veis the equilibrium volume of the basin.For a one-element model for the back-barrier basin,the total water volume of the basin below MHW,V,which is the sum of the channel volume(below LHW)Vcand the tidal prism P,can be expressed as V=Vc+P.

A.3.Sediment transport and morphological change

As explained in the previous section,the sediment transport processes are described by the mass balance equation for sediment for each element aggregated in time(Eq.(7)).Following the schematisation shown in Fig.1,these equations are

where J is the sediment transport,with the subscripts indicating the transport directions,in which the subscripts f,c,d,and o mean the intertidal flat,channel,ebb-tidal delta,and outside area,respectively.Jcfis thus the sediment transport from the channel to the intertidal flat.FBf,FBc,and FBdare sediment exchange flux between the bed and the water column of the intertidal flat,channel,and ebb-tidal delta,respectively.Because of the aggregation in time,the transports are governed by the diffusion/dispersion with the tidal flow as mixing agent:

where c is the suspended sediment concentration,and its subscript indicates the morphological element.δ is the horizontal exchange coefficient,which depends on the dispersion coefficient,the area of the cross-section connecting the two elements,and the distance between the two elements(see Wang et al.,2008).Its subscripts have similar meanings as those of J.

FBis the sediment exchange flux between the bed and water column,representing the erosion minus sedimentation.FBf,FBc,and FBdare formulated as follows:

where Af,Ac,and Adare the horizontal areas of the intertidalflat,channel,and ebb-tidal delta,respectively;cfe,cce,and cdeare the equilibrium sediment concentrations of the intertidal flat,channel,and ebb-tidal delta,respectively;and wf,wc,and wdare the vertical sediment exchange velocities of the intertidal flat,channel,and ebb-tidal delta,respectively,which are proportional to the settling velocity of sediment.Note that the vertical exchange velocity is not necessarily the same for all three elements even if only a single fraction of sediment is present.As discussed earlier,the equilibrium sediment concentrations for the elements are calculated by comparing the morphological state with its equilibrium state:

where Cf,Cc,and Cdare coefficients for different elements,with the dimension of sediment concentration;and n is a power similar to that in a power-law sediment transport formula(Wang et al.,2008),with the subscripts indicating the morphological elements.Note that there is a difference in the formulations between the elements of the intertidal flat and ebb-tidal delta with a sediment volume as a morphological state variable and the channel element with a water volume as a state variable(sediment volume relationships are the inverse of the water volume relationships).

The sediment exchange between the bed and water column(Eqs.(A.13)-(A.15))also determines the change of the morphological state variable of the corresponding element,with dVf/dt=-FBf,dVc/dt=FBc,and dVd/dt=-FBd.If all the morphological elements in the whole system are in equilibrium,FBshould vanish for all elements.Thus,the sediment concentration is equal to its equilibrium value for each of the elements,with cf=cfe,cc=cce,and cd=cde.Eqs.(A.7)-(A.9)yield Jcf=Jdc=Jod=0,i.e.,all residual sediment transports vanish,followed by cf=cc=cd=cofrom Eqs.(A.10)-(A.12),and also cfe=cce=cde.According to Eqs.(A.16)-(A.18),cfe=Cf,cce=Cc,and cde=Cdif all the elements are in equilibrium(their volumes are equal to their equilibrium values),and then,Cf=Cc=Cd=co.Thus,the coefficients in Eqs.(A.16)-(A.18)should be the same and equal to the sediment concentration at the open boundary(the outside area).This constant is called the global equilibrium concentration cE,with Cf=Cc=Cd=co=cE.Any deviation from this rule implies a redefinition of the morphological equilibrium state.