This study tests the hypothesis that horizontal density gradients have the potential to significantly contribute to the accumulation of suspended particulate matter (SPM) in the Wadden Sea. It is shown by means of long-term observations at various positions in the Wadden Sea of the German Bight that the water in the inner regions of the Wadden Sea is typically about 0.5–1.0 kg m−3 less dense than the North Sea water. During winter this occurs mostly because of freshwater runoff and net precipitation; during summer it occurs mostly because of differential heating. It is demonstrated with idealized one-dimensional water column model simulations that the interaction of such small horizontal density gradients with tidal currents generates net onshore SPM fluxes. Major mechanisms for this are tidal straining, estuarine circulation, and tidal mixing asymmetries. Three-dimensional model simulations in a semienclosed Wadden Sea embayment with periodic tidal forcing show that SPM with sufficiently high settling velocity (ws = 10−3 m s−1) is accumulating in the Wadden Sea Bight because of density gradients. This is proven through a comparative model simulation in which the dynamic effects of the density gradients are switched off, with the consequence of no SPM accumulation. These numerical model results motivate future targeted field studies in different Wadden Sea regions with the aim to further support the hypothesis.
a. Wadden Sea as suspended matter sink
Interest in suspended particulate matter (SPM) has increased during the last decades not only because contaminants like lead and organic pollutants [e.g., polychlorinated biphenyls (PCBs) or tributyltin oxide (TBT)] bind to SPM rather than dissolve in the water phase, but also because they absorb light (and thus hinder primary production). In areas of high SPM concentration like the Wadden Sea, pelagic primary production decreases with respect to the open sea due to high water turbidity (Van Beusekom et al. 2001).
For this reason, the SPM concentration in the Wadden Sea, its variability, and its transport has been investigated in the past 50 yr with many experimental studies and modeling exercises. One of the most interesting features is that any SPM satellite scene (see, e.g., Fig. 1) shows a high SPM concentration inside and a low concentration outside the Wadden Sea. This well-known feature alone should give rise to an outward-directed SPM flow due to turbulent dispersion. However, the Wadden Sea is more or less stable on a time scale of decades or even centuries. Measurements, on the other hand (Pejrup 1988a, b; Eisma 1993; Puls et al. 1997; Townend and Whitehead 2003), often show an inward transport (resulting in a sea bottom rise of 0.07 to 8 mm yr−1) that may be altered by human influences like the construction of dikes, land reclamation, or coastal lagoon closure (Afsluitdijk in the Netherlands separating the Lake Ijssel from the North Sea, and the Hindenburg and Rømø dams closing the Sylt–Rømø Bight to the south and north, respectively; see Fig. 2).
Observational evidence shows that suspended matter deposition is inhomogeneous in time and in space inside the Wadden Sea. Suspended matter is accumulated in summer but eroded during winter or during severe storms (e.g., Lumborg and Pejrup 2005). The deposition rate may be high on tidal flats but low or negative in the tidal channels (Dyer 1994; Pedersen and Bartholdy 2006). Thus, averages like annual accumulation rates are difficult to estimate. Moreover, the annual sediment budget can be largely determined by extreme storm events (see Lumborg and Pejrup 2005). A very effective method therefore derives deposition rates from sediment cores by measuring activities of radioactive elements. However, not many areas are suited for this purpose because the Wadden Sea sediment is often reworked in the upper 10 cm by benthic macrofauna (Andersen et al. 2000).
Eisma (1993) estimated that the Wadden Sea sediment growth, in general, keeps pace with sea level rise and land subsidence; he reports sedimentation rates of 10–28 mm yr−1 for the Dollart, Jade Bay, and Leybucht in the German Bight. De Haas and Eisma (1993) reported 1–2 mm yr−1 (with variations up to 8 mm yr−1) in the Dollart. Postma (1981) estimates that the tide entering the Wadden Sea leaves 3.5 × 106 t yr−1 of suspended matter behind, which amounts to an annual average growth of 0.27 mm (with 8500-km2 area for the Wadden Sea and 1500 kg m−3 density for the SPM). Pejrup (1997) measured with 210Pb dating an accumulation rate in the Sylt–Rømø Bight of 58 000 t yr−1, 63% of which are due to the water exchange with the North Sea. With an area of 400 km2 for the bight, this amounts to only 0.06 mm yr−1. Andersen et al. (2000) and Andersen and Pejrup (2001) measured rates of 5–12 mm yr−1 for a specific site in the Sylt–Rømø Bight demonstrating that the accretion rates can differ by two orders of magnitude inside the same basin.
However, there is also some evidence of sediment export: Flemming and Nyandwi (1994) argue that the combined influence of sea level rise and the construction of dikes leads, by an increase in the mean tidal velocity, to a net export of fine-grained material. Thus one or more mechanisms must counterbalance or even exceed the outward-directed turbulent dispersion of SPM.
For such a mechanism, several suggestions have been made:
Settling lag (Postma 1954; Van Straaten and Kuenen 1958; Postma 1961; Bartholdy 2000): Sedimented material is taken into suspension by the flood current at a certain critical water velocity and transported onshore. Once the velocity falls below the critical value (because of changing tide or decreasing current velocity toward the mainland), the particles do not settle immediately, but take some time to reach the ground, moving farther inward than a symmetric water velocity pattern would allow. The basic assumption of the settling lag mechanism (and also the second mechanism) is a nonzero sinking velocity for SPM and an inwardly decreasing water velocity that introduces the asymmetry.
Scour lag (Postma 1954; Van Straaten and Kuenen 1958; Postma 1961; Bartholdy 2000): During ebb, a higher current speed is required to erode the particles from the bottom than the current speed at which the particles have settled during the flood before. This is due to the fact that critical shear stresses above which sedimentation ceases are typically lower than critical shear stresses above which erosion takes place. Bartholdy (2000) has shown for the northern Danish Wadden Sea that the scour lag is much more important than the settling lag.
Asymmetric tidal water level curve and current velocities (Groen 1967; Dronkers 1986a, b): Lumborg and Windelin 2003 explain in a qualitative way the effects of asymmetric current profiles on the sediment transport, and two asymmetries are taken into account. First, the change in current velocity during high water is much slower than that during low water. Therefore, the suspended sediment has more time to settle down on the tidal flats than at low tide where water is mainly in the deeper tidal channels. Second, the current velocity is higher during flood than during ebb; thus the sediment is carried farther inward. Even if the tidal wave has a symmetric shape outside the Wadden Sea, it may be distorted when entering the tidal channels. Lumborg and Windelin (2003) argue that the tidal wave once entering one of the Wadden Sea Bights is distorted because of decreasing water depths, resulting in higher current speed during flood than during ebb. Each of these effects may lead to inward transport of suspended matter. However, Ridderinkhof (2000) found that, in the Ems-Dollard estuary, the velocity profiles are opposite to the ones reported from the western Dutch and northern Wadden Seas; that is, the ebb current velocities are stronger than the flood velocities. Moreover, the change in current velocity is faster during flood tide than during ebb. Thus, the inward transport mechanism caused by these kinds of asymmetry is not likely to play a major role.
Flocculation (Pejrup 1988a; Dyer 1994): Higher suspended sediment concentration in the inner parts of the Wadden Sea leads to coagulation into sediment flocs with higher settling velocity. Thus, these particles settle faster than those near the tidal inlets with a lower SPM concentration, yielding a net inward transport. In this mechanism, the basic assumption is an increasing suspended matter concentration toward the mainland, which may have been established by one of the other mechanisms discussed here. Thus flocculation may have the potential to amplify other mechanisms.
Stokes drift (Dyer 1988, 1994; Stanev et al. 2007), a second-order effect of the tidal wave: The tidal wave is distorted because of decreasing water depth on its way toward the inner part of the Wadden Sea. The Stokes drift tends to carry water into the tidal-wave-traveling direction; in order to maintain continuity, horizontal residual circulation cells may develop. Dyer (1994) reports that Stokes drift can cause inward as well as outward net transport in estuaries depending on the relationship between the cross-sectional bathymetry and the tidal water level elevation.
No quantitative study has hitherto established a ranking of all of these effects in the Wadden Sea. Furthermore, none of these net SPM transport mechanisms has considered horizontal density differences in the Wadden Sea as potential driving force for suspended matter transport.
In this paper, we investigate the hypothesis that horizontal density gradients can contribute to net SPM accumulation in the Wadden Sea. First, we describe monitoring pile stations in the Wadden Sea (section 2a) and discuss time series observations of salinity, temperature, and density recorded by instruments mounted on these piles (section 2b). Afterward, density-gradient-related mechanisms potentially driving net SPM transport are discussed in section 2c. Numerical modeling tools used for a quantitative analysis are introduced in section 3, including a turbulence closure model (section 3a) and a coastal dynamics model (section 3b). One-dimensional model experiments with a prescribed horizontal salinity gradient are then used for demonstrating the mechanisms leading to net SPM transports down the salinity gradient (section 4). Three-dimensional model simulations of a semienclosed tidal embayment with freshwater runoff are then used to show the effect of density gradients on the SPM accumulation in the Wadden Sea (section 5). A discussion of the results including suggestions for further testing of the hypothesis is presented in section 6.
2. Density gradients in the Wadden Sea
a. Materials and methods
Time series of hydrographic parameters were measured by GKSS (i) at two locations in the northern (Rømø Dyb) and the southern (Lister Ley) tidal channels of the Sylt–Rømø Bight during 1996–97, (ii) at a location on the outer reach of the Accumersieler Balje, the tidal inlet between the East Frisian islands Langeoog and Baltrum, in the months of May to November since 2000, and (iii) at a location in the inner parts of the Hörnum Deep between the islands Sylt and Föhr in the North Frisian Wadden Sea since 2002 (see Fig. 2). The setup of these measuring stations was as follows: An iron pile with a diameter of 0.40 m was fastened several meters into the ground close to the 3-m low-water line. Underwater probes for conductivity, water temperature, optical transmission, current velocity, gauge level, and wave height as well as a meteorological station measuring air temperature, air and pressure, wind speed and direction, global irradiation, and precipitation were mounted onto these piles. The underwater sensors were installed 1 to 1.5 m above the seafloor. The instrument specifications are given in Table 1. Power supply was granted by two solar modules. Ten-minute averages of 2-Hz time series were collected at the pile station and transmitted to a mainland station where the data were stored. Every hour they were further transmitted to GKSS via phone connection. Salinity was calculated from conductivity and temperature using the United Nations Educational, Scientific and Cultural Organization (UNESCO) formula of Rathlev (1980). Data from all four GKSS pile stations are available online at http://tsdata.gkss.de.
Different from the GKSS piles, the Spiekeroog pile (built by the University of Oldenburg) is a permanent time series station that was set up in summer 2002 at the tidal inlet between the two East Frisian islands Langeoog and Spiekeroog (see Fig. 2) at an average water depth of 13 m. Its construction withstands ice thicknesses of up to 50 cm. The mechanical structure consists of a pile having 35.5-m length and a diameter of 1.6 m, driven to 10-m depth into the sediment. The station provides continuous data on temperature and conductivity all year-round taken from five different depths (1.5, 3.5, 5.5, 7.5, and 9 m above the seafloor) throughout the water column as well as pressure at 1.5 m above the seafloor. The instrument specifications are given in Table 1. The sensors in the pile are located in specific tubes with a horizontal orientation along the main tidal current direction, thus providing a continuous flow of seawater along the sensors. Meteorological data include temperature, pressure, wind speed and direction, air humidity, downwelling irradiance, as well as sky and water leaving irradiance. In addition, spectral light transmission, nutrients, oxygen, yellow substance, and methane are measured below the sea surface. A bottom-mounted ADCP is located close to the station. Power supply is granted by a wind energy converter and solar panels. A 10-min average of the data is sent to a receiving station on Spiekeroog island where the data are fed via telephone network to the University of Oldenburg five times a day.
For quality assurance, depth profiling measurements using precise CTD probes are either carried out on the station or aboard an anchored research vessel. The sensors are calibrated by the manufacturer and checked in the laboratory directly before deployment at the time series station. Salinity is calculated based on the equation of the state (Fofonoff and Millard 1983). Details on the time series station can be found in Reuter et al. (2007, manuscript submitted to Meas. Sci. Technol.). The real-time data are available on the Internet at http://las.physik.uni-oldenburg.de/wattstation.
The data coverage of all five pile stations is shown in Fig. 3. For the Langeoog and Hörnum Deep piles, the salinity and temperature data have been averaged over one month using all available data from 2000 to 2006. The error bars in Fig. 4 have been calculated using the bootstrapping method.
b. Long-term observations in the Wadden Sea
The idea that the Wadden Sea water may be lighter than adjacent North Sea water arises from salinity and temperature observations, salinity having a much larger impact on density than temperature. At 10°C and 31-psu salinity, an increase of 2°C results in a density decrease by 0.35 kg m−3, whereas a decrease of 2 psu leads to a density decrease of 1.6 kg m−3. On the one hand, precipitation has a larger impact in the Wadden Sea than in other areas because of the region’s shallowness; for example, a rainfall event of 20 mm adds, in an area with only 2-m average water depth, as much as 1% to the total water mass and decreases salinity accordingly. This effect is intensified by the fact that the water depth is not uniform in the Wadden Sea. When precipitation occurs at low tide, fresh rainwater falls on the dry tidal flats and produces a water body of very low salinity. Because of the strong currents, vertical mixing is very effective and results in a horizontal salinity gradient from the tidal inlets toward the shore. This gradient is subsequently reduced by shear dispersion between the Wadden Sea water body and the North Sea within a few days.
Additionally, the Wadden Sea receives freshwater from sluice gates to the mainland draining the low-lying areas behind the dikes. Thus, even if there is no permanent freshwater source, the Wadden Sea water usually has a lower salinity than the adjoining North Sea water body. This is demonstrated in the upper panel of Fig. 4. In the East Frisian Wadden Sea, salinity differences between high and low water drop from >1.6 psu in winter and spring months to <0.5 psu in summer; in the Hörnum Deep, the change is from approximately 0.5 to <0.1 psu. The Spiekeroog data show the entire seasonal cycle and vary from 0.8 psu in January to about 0 in August. Although these two piles are situated in a very similar environment, the salinity differences of the Langeoog pile are much larger, which is probably due to the difference in its location. The data from the Sylt–Rømø Bight fall in between the range defined by Langeoog and Spiekeroog values with a similar seasonal behavior. The seasonal variation reflects the variability in rainfall and evaporation.
Using both freshwater inflow from precipitation and sluices and the salinity difference, mean water exchange time scales can be calculated. Assuming an annual rainfall of 800 mm (as is average for northern Germany), an area of 90 km2 for the Wadden Sea Bight between Langeoog, Baltrum, and the mainland, an area of 140 and 75 km2 for the catchment areas of the sluices Accumersiel and Bensersiel, respectively, and further assuming that 50% of the precipitation on land is lost by evaporation, we estimate an average freshwater inflow of 433 000 m3 day−1. This amounts to 2.76 ppt of the water volume at the mean tidal level, thus reducing the mean salinity of roughly 31 psu by about 0.085 psu day−1. If, furthermore, we consider an average salinity difference of 0.4 psu between the North Sea and the Wadden Sea waters, we estimate that the mean adjustment time is about 4.7 days. The average flushing time, that is, the ratio of the tidal prism to the inflow is, however, much longer and amounts to about several months. Of course, these time scales are strongly modified by storms as well as by the spring–neap cycle, but they are valuable as rough estimates.
The second process that alters the Wadden Sea water density is temperature. As can be seen in the middle panel of Fig. 4, Wadden Sea water is warmer than North Sea water from March to August (with a maximum difference of about 2 K). The reason for the temperature difference oscillation may be the varying day length: from March to September, the days are longer than the nights, so that the warming of shallow water and tidal flats during daytime dominates over the cooling during the night. From the end of September to March, the opposite effect (longer nights) leads to a stronger cooling in the Wadden Sea and thus to an inverse temperature gradient.
Both influences on the water density can be combined to calculate mean density differences between high water and low water (lower panel of Fig. 4). All pile measurements demonstrate that the Wadden Sea water is usually lighter than the North Sea water by about 0.4 kg m−3. The pile at Langeoog shows this influence more clearly than the other ones because of the position of the pile, which is ideal for the detection of Wadden Sea as well as North Sea waters. This general property of an onshore to offshore density gradient may give rise to additional tidal asymmetries, which must be taken into account when considering the mechanisms carrying suspended matter into the Wadden Sea. The potential impact of this is discussed in the next section.
c. Potential dynamic impact of density gradients
The processes discussed in the previous section have the consequence that coastal waters in regions of freshwater influence (ROFIs; see e.g., Sharples and Simpson 1995) are substantially less dense than the offshore waters. Where tides are present in such regions, coastal transports may be dominated by a nonlinear interaction between the tidal currents and the horizontal density gradients. Simpson et al. (1990) were the first to systematically investigate these effects. During ebb tide less dense coastal waters are sheared over denser offshore waters, creating a tendency to statically stabilize the water column. Depending on the ratio between the stratifying shear and the shear-generated turbulence, the water column may become stratified or not. This interaction of forces is called tidal straining (see Simpson et al. 1990). During flood, the opposite happens: denser offshore waters are sheared over less dense onshore water, increasing the potential energy of the water column, which is then released into turbulent kinetic energy (TKE), subsequently leading to enhanced vertical mixing. Thus, during flood the water column has the tendency to become more mixed. When a certain balance between the tidal force and the horizontal pressure gradient is given, the water column may be stratified after ebb and well-mixed after flood. This so-called strain-induced periodic stratification (SIPS) has been intensively studied in Liverpool Bay by Simpson et al. (1990), Sharples and Simpson (1995), Rippeth et al. (2001), and Simpson et al. (2002).
In tidal estuaries where horizontal salinity gradients and tides interact, processes similar to SIPS are present. Although it may be that an estuary is permanently stratified, SIPS may act near the bed, leading to near-bed destabilization during flood and to near-bed stabilization during ebb. As Jay and Musiak (1994) pointed out, the systematic change of increased turbulence during flood and suppressed turbulence during ebb leads to a tidal asymmetry also in eddy viscosity. Stronger mixing during flood transports more momentum toward the bed, leading to a transport asymmetry that results in a near-bed residual upstream transport. Since SPM concentrations increase toward the bed because of the interaction between settling and vertical mixing, a net upstream transport of suspended matter typically results. The upstream end of this transport is given by the limit of the salt intrusion, such that suspended matter is accumulating there as an often clearly visible estuarine turbidity maximum (ETM). The ETM formation is actually supported by a number of other processes, such as, for example, the estuarine circulation (dense water tongue creeping upstream against the mean flow; see Postma and Kalle 1955; Festa and Hansen 1978; Dyer 1988), or the fact that suspended matter is mixed higher up into the water column during flood, due to increased vertical mixing (see Geyer 1993). Burchard and Baumert (1998) have compared the efficiencies of these three mechanisms with the aid of an idealized numerical model study and found that the tidal asymmetry suggested by Jay and Musiak (1994) was the most efficient. It should be noted that a number of other ETM mechanisms not related to density differences have been identified, such as barotropic tidal asymmetries (see, e.g., Allen et al. 1980) or topographic features such as bends, headlands, or sills (see, e.g., Fischer et al. 1979), leading at times to multiple ETMs (see Rolinski 1999).
Postma (1954) and Zimmerman (1976) discuss in detail vertical and horizontal salinity distributions in the Marsdiep, a tidal channel connecting the western Dutch Wadden Sea with the southern North Sea. Because of substantial freshwater supply from Lake Ijssel, high to low water salinity differences are typically of the order of 2 psu, and thus higher than in the East and North Frisian Wadden Sea (see Fig. 4). Since the vertical salinity profiles are well mixed, Postma (1954) concluded that gravitational circulation does not play a significant role in the Marsdiep. The other density-related potential driving mechanisms for residual SPM transport as described by Geyer (1993) and Jay and Musiak (1994) have not yet been discussed for the Wadden Sea.
The question that is investigated in the present work is whether the superposition of strong tides with small horizontal density gradients, as is typical for the Wadden Sea (see section 2), has the potential to generate significant net transports of suspended matter by means of the density-related mechanisms discussed above.
3. Modeling tools
a. Turbulence closure model
During the last decades, two-equation turbulence closure models have proven to be most suitable for reproducing the interaction between turbulence and mean flow in buoyancy-driven coastal flow dynamics. These statistical models, which use one dynamic equation for the TKE and one other related to the integral length scale of turbulence, are based on statistical treatment of the Navier–Stokes equations, and do thus contain a number of empirical parameters. Those parameters have physical meanings and are thus adjustable by means of laboratory and model [direct numerical simulation (DNS)] experiments. For more details of the parameter determination, see Umlauf and Burchard (2003). In the closure approach we apply here, second-moment dynamics are reduced to algebraic relations, a simplification that makes such closures suitable for implementation into three-dimensional models; for details, see Burchard and Bolding (2001). In several applications, these models have proven to agree quantitatively well with turbulence observations in coastal and limnic waters; see, for example, the investigations by Burchard et al. (2002), Simpson et al. (2002), and Stips et al. (2002, 2005.
The closure model that is applied in the present study is the k–ɛ model with transport equations for the TKE, k, and the turbulence dissipation rate, ɛ. From k and ɛ, the integral length scale can be calculated by assuming equilibrium turbulence spectra. As second-moment closure, we apply the model suggested by Cheng et al. (2002).
b. Coastal dynamics model
The dynamics of bathymetry-guided near-bed processes such as they occur in the Wadden Sea are best reproduced by numerical models with bottom-following coordinates. The major advantage of such model architecture as compared to geopotential coordinates with steplike bottom approximation is that for bottom-following flows, the advective fluxes across vertical coordinates, and thus the associated discretization errors, are minimized (see Ezer and Mellor 2004; Ezer 2005). Furthermore, bottom-following coordinates allow high near-bed resolution also in domains with large depth variations (see Umlauf and Lemmin 2005). Compared to these advantages, the pressure gradient discretization error due to isopycnals crossing sloping coordinates (see, e.g., Haney 1991; Kliem and Pietrzak 1999), an error typically associated with bottom-following coordinates, is of minor importance for the energetic Wadden Sea.
The three-dimensional General Estuarine Transport Model (GETM; see Burchard and Bolding 2002; Burchard et al. 2004), which has been used for the present numerical study, combines the advantages of bottom-following coordinates with the turbulence module of GOTM (see section 3a). A numerical feature that supports high numerical resolution computations are high-order advection schemes, as described by Pietrzak (1998). The basic concept behind these schemes is to apply well-tested monotone [total variation diminishing (TVD)] one-dimensional methods in a directional split mode. The scheme used here is a third-order monotone ULTIMATE QUICKEST method; see Leonard (1991) and Pietrzak (1998) and for numerical tests as well as Burchard and Bolding (2002). The major advantage of these schemes is that they retain frontal features without creating spurious oscillations and negative concentrations. Furthermore, high spatial resolution is obtained by the possibility of running GETM on parallel computers.
Simulations of scenarios with occurrence of intertidal flats are supported because of a dynamic drying and flooding algorithm. In contrast to many other models, GETM is not explicitly shutting off transports out of grid boxes in which the water depth falls below a critical depth. Instead, the dynamics are retained in such a way that it gradually reduces to a balance between external pressure gradient and friction. For details, see Burchard et al. (2004).
GETM has been successfully applied to several coastal, shelf sea, and limnic scenarios [see, e.g., Stanev et al. (2003a, b) for turbulent flows in the Wadden Sea, Stips et al. (2004) for dynamics in the North Sea, Burchard et al. (2005, 2007) for studies of dense bottom currents in the western Baltic Sea, Banas and Hickey (2005) for estimating exchange and residence times in Willapa Bay in Washington State, and Umlauf and Lemmin (2005) for a basin-exchange study in the Lake of Geneva].
4. One-dimensional analysis
a. One-dimensional dynamic equations
Let us consider a water column of depth H with the upward-directed coordinate z, representing the vertical structure of a water body of horizontally infinite size. An external pressure gradient force is oscillating periodically with period T, superimposed on a constant horizontal density gradient ∂xρ aligned with the direction of the external pressure gradient. With the buoyancy
and neglecting earth’s rotation, the momentum equation is then of the form
where the first term on the right-hand side represents the effect of the internal pressure gradient, and pg is a nondimensional periodic external pressure gradient function with period T chosen such that
with the vertical mean velocity amplitude umax. This construction guarantees that the vertically and tidally integrated transport ∫T0 U(t) dt is zero (see Burchard 1999 for details). In (1), g denotes the gravitational acceleration and ρ0 the constant reference density, and in (2), υt denotes the eddy diffusivity. The effect of earth’s rotation can be neglected for this idealized one-dimensional study, since the external Rossby radius is much larger than the dimensions of the study area and stratification is weak.
The transport equation for salinity is given as
where the term on the right-hand side is a relaxation term. Without this term, salinity would run away because of residual horizontal buoyancy transport. In real estuaries, this relaxation to a certain tidal mean buoyancy is given by complex three-dimensional mixing processes that are missing here. Thus the relaxation term can be considered as a parameterization for these lateral processes. In (4), υ′t denotes the eddy diffusivity. The density ρ (and thus the buoyancy b) is then calculated by a linear equation of state; see Eq. (13). The horizontal salinity gradient is prescribed such that the buoyancy gradient ∂xb (see above) can be calculated.
A simple but general suspended matter equation is given by
with the nondimensional suspended matter concentration C and the positive constant settling velocity ws. There are no fluxes of suspended matter through the bottom and the surface, that is, at low near-bed turbulence; suspended matter is accumulating near the bed as a kind of fluff layer and is mixed up into the water column when the near-bed turbulence is increasing again. Furthermore, only low-concentration, noncohesive sediments with no feedback to the buoyancy of the flow are considered here, such that a nondimensional view with constant settling velocity is fully sufficient here.
As turbulence model, a two-equation k-ɛ model is chosen here. The TKE equation is here of the form
with the TKE k, the dimensionless turbulent Schmidt number σk, and the turbulent dissipation rate ɛ, for which an empirical equation of a form similar to (6) is used (see, e.g., Burchard 2002). The terms on the right-hand side of (6) represent shear production, buoyancy production, and dissipation.
The three potential, density gradient related, driving mechanisms for residual SPM transport discussed in section 2c are included in this one-dimensional model. Gravitational circulation is included by the combination of the internal and external pressure gradient on the right-hand side of (2). SIPS is considered because of the combination of salinity advection in (4) and the stratification effects included in the turbulence model [e.g., the buoyancy production term in (6)] and thus it is included in the eddy viscosity and diffusivity in (2) and (4). Finally, the tidal mixing asymmetry (Geyer 1993) is included in (5) via the stratification dependence of the eddy coefficient for SPM.
It should be noted that in this idealized case, the system of hydrodynamic and turbulence equations depends on only two nondimensional parameters, which are the inverse Strouhal number (see Baumert and Radach 1992),
and the horizontal Richardson number (see Monismith et al. 1996),
An increase in depth or a decrease of tidal velocity amplitude would specifically lead to an increased Rx, and with this to a more pronounced tidal asymmetry.
b. One-dimensional process study
For a one-dimensional model study, we choose specifications that are typical for tidal channels in the Wadden Sea: a water depth of H = 20 m, a vertical mean velocity amplitude of umax = 1 m s−1, a tidal period of T = 44 714 s (in accordance with the M2 constituent), and a horizontal salinity gradient of ∂xS = −1.4 × 10−4 psu m−1 such that the x coordinate increases in onshore direction. Salinity is relaxed to 30 psu, with a relaxation constant equal to the tidal period. Sea level variations are not considered here. The simulation is initialized from some neutral conditions (such as zero velocity, minimum turbulence, constant salinity) and is carried out for 10 tidal periods, of which the last 5 ones turned out to be in periodically steady state. Results are shown for the last tidal period only. The high tide to low tide differences in salinity (see Fig. 5) are about 2 psu and thus do agree well with typical winter Wadden Sea conditions (see Fig. 4).
The SIPS mechanism is clearly seen in the right panel of Fig. 5: at full flood, the surface salinity exceeds the bottom salinity because of stronger positive salt advection near the surface. This is an unstable state, which lasts until the end of flood. An interesting feature during ebb is the stratification maximum during early ebb (t/T = 0.65), which is explained by the time lag between the stratifying shear and the destratifying shear-generated vertical mixing. This has the consequence that the eddy viscosity is substantially smaller during full ebb than during full flood (see the left panel in Fig. 6 where we show profiles shortly after full flood and full ebb) in order to allow the flow to adjust to the flow conditions. During flood the flow velocity is thus more homogeneous in the bulk of the water column, with higher speeds near the bottom and lower speeds near the surface, when compared to the ebb flow (see Fig. 6, middle panel). This has the consequence that near-bed tidal mean flows are directed toward the coast (see Fig. 6, right panel).
For SPM, two different settling velocities are considered, with values of ws = 10−3 m s−1 (fast) and ws = 10−4 m s−1 (slow). The left panels of Figs. 7 and 8 show SPM profiles shortly after full flood and shortly after full ebb. One can clearly see that during flood (because of high turbulence) the profiles are more homogeneous than during ebb (see left panels of Figs. 7 and 8).
The tidally averaged SPM transport can be decomposed into two contributions. By introducing the tidal averaging operator
and denoting the averages and fluctuations of velocity and SPM profiles,
Tidal mean SPM transport profiles 〈uC〉 (see the right panels of Figs. 7 and 8) are directed landward near the bed and seaward near the surface, with a landward-oriented residual. For the fast-settling SPM, about half of this transport is due to the overlay of tidally averaged velocity and SPM profiles uC, and the other 50% are due to an on-average positive correlation of ebb currents with low SPM concentrations and flood currents with high SPM concentrations. For the slow-settling SPM, most of the tidally averaged SPM transport is explained by an overlay of averaged velocity and SPM profiles, except near the surface.
The tidally and vertically integrated suspended matter transport,
is Cint = 6235 m2 for the fast-sinking SPM and Cint = 945 m2 for the slowly sinking SPM. Although these are still small fractions of the total suspended matter transport, which are moved onshore and offshore during each flood and ebb tide (about 2.2% for the fast and 0.33% for the slow SPM), this shows that there is a buoyancy-related net transport mechanism for a wide range of SPM classes.
As discussed in section 2c, these processes are well known for estuaries and river plumes. However, for Wadden Sea dynamics they have so far not been considered. As shown in this section, they may nonetheless considerably contribute to a net onshore transport of SPM also in Wadden Sea regions.
5. Three-dimensional model studies
The tidal asymmetries in stratification and current shear are so small that it is difficult to detect them by means of field observations. In order to study the potential effect of horizontal density differences in a real tidal basin, we therefore apply a three-dimensional numerical model for an idealized situation to the Sylt–Rømø Bight in the Wadden Sea area in the vicinity of the Danish–German border. We chose this specific area because artificial dams in the south and north of the domain create a semienclosed tidal embayment with only one entrance (see the right inset in Fig. 2). Because of these idealized conditions, many numerical model studies have been carried out in this domain (see, e.g., Dick 1987; Burchard 1998; Schneggenburger et al. 2000; Lumborg and Windelin 2003; Lumborg and Pejrup 2005; Villarreal et al. 2005).
a. Model setup
For the three-dimensional simulations in the Sylt–Rømø Bight, GETM (see section 3b) is used in baroclinic mode. The model is forced through prescribed surface elevations at the open boundaries in the north, west, and south (see Fig. 9). These surface elevations have been composed on the basis of harmonic analysis of the dominant M2 tide as computed by the operational model of the German Federal Maritime and Hydrographic Agency (BSH) with a period of T = 44 714 s. Other tidal constituents have not been considered in order to obtain a period model simulation (Burchard 1995). Freshwater flux from river runoff is specified at two locations at the eastern and northeastern end of the Sylt–Rømø Bight at constant rates of 10 m3 s−1 for each of the two small rivers Brede Å and Vida Å (see Fig. 9). Together with a constant precipitation rate of P = 0.1 m month−1 this results in a total freshwater flux of 36 m3 s−1 into the inner Sylt–Rømø Bight. These freshwater flux values are typical for a winter month. The inner Sylt–Rømø Bight is defined here as the area composed of budget boxes I (south) and II (north) as indicated in Fig. 9. Heat fluxes are not considered here, such that the simulation is carried out with constant temperature. Surface momentum fluxes are also set to zero, such that wind and surface wave effects are not considered in this idealized simulation.
The model is initialized with zero surface elevation η (however, η = −H + hmin for intertidal flats, where H is the depth below zero and hmin = 0.05 m is the minimum water depth), zero velocities, and minimum values for the turbulence quantities. Salinity is initialized with 28 psu at the mouth of the Vida Å, increases linearly with distance from there to 32 psu at a distance of 26 km, and has a constant value of 32 psu at distances larger than 26 km from the mouth of the Vida Å. SPM is initialized to unity in the entire domain after 45 tidal periods of model spinup. River runoff has zero salinity and SPM concentration as well as the precipitation.
At the open boundaries, salinity is relaxed to 33 psu during inflow with a relaxation time scale of τ = Δx/un, where Δx = Δy = 200 m is the horizontal resolution of the model grid and un is the vertically averaged normal component of the velocity at the boundary. At the open boundaries, SPM is relaxed to unity with the same time scale. Using fixed boundary values for salt and SPM also under inflow conditions instead (equivalent to zero-relaxation time scale) would lead to the propagation of unrealistically sharp fronts into the domain.
The water column is resolved with 20 equidistant σ layers, such that the vertical model resolution varies between about 1.5 m at the deepest point (29.1 m) and 0.0025 m over tidal flats. Each tidal period is resolved with 104 barotropic and 103 baroclinic time steps. Advection terms for salt, momentum, and SPM are discretized with the positive-definite ULTIMATE QUICKEST scheme (see section 3b). For the vertical advection, a time-split scheme is additionally used in order to guarantee numerical stability. The settling velocity of SPM is linearly reduced to zero between the critical water depth hcrit = 0.2 m, and the minimum water depth hmin = 0.05 m, in order to avoid too many iterations of the vertical advection scheme in tidal flat regions.
Two types of simulations have been carried out, one with consideration of density differences and one without. Potential density ρ is calculated from salinity S by means of the following linear equation of state:
where the haline expansion coefficient β = ∂Sρ has the following values:
b. Simulation results
Both simulations with and without consideration of density gradients were first spun up from initial conditions for 45 tidal periods without calculation of SPM. Between the two budget boxes in the inner bight, a west–east section is defined, starting from the northwestern tip of the island of Sylt and cutting through the southern branch of the tidal channel (see Fig. 9). At the deepest point of this section with a water depth of 25 m, position P, time series data are extracted.
It can be seen from Fig. 10, which shows time series of surface and bottom salinity at position P, that during flood unstable stratification is dominant, and during ebb stable stratification is dominant. The salinity amplitude at position P is about 2 psu, with values oscillating around 30 psu. This is similar to the results of the one-dimensional simulations (see Fig. 5). The time series for salt at position P also shows that a quasi-periodic state was reached after 55 tidal periods, with a small trend in salinity remaining. The latter can be explained by a water exchange time of the inner Sylt–Rømø Bight of about 265 days, which results from division of its average volume by the freshwater fluxes and which is much longer than the duration of 45 tidal periods (about 23 days). However, the density gradients are basically in periodic state after 45 tidal periods, as can be seen from Fig. 10.
Strong horizontal salinity gradients are built up in the regions of the river runoff locations (see Fig. 11), specifically in the fairly narrow northern Sylt–Rømø Bight. This is supported by the net precipitation, which also generates very low salinities on the sandbanks southwest of the island of Rømø.
The significant differences between ebb and flood currents are shown in Fig. 12 for the area northeast of the island of Sylt. The instants of flood and ebb are defined such that for ebb flow the instant of maximum current speed along the west–east section (see Fig. 9) is chosen, and for flood an instant shortly after maximum flood is chosen that has a similar maximum velocity as the maximum ebb along this section. This is done in order to provide better intercomparison of ebb and flood results for this section. The curved flows around the northeastern end of Sylt are subject to inertial effects, such that flood currents generate a lee zone south of the northeastern end and ebb currents generate a lee zone north of it. These tidal asymmetries are more clearly visible in the tidal residual velocities, defined as
with the vertically integrated transport vector (U, V) and the water depth D (see Fig. 13).
This can also be clearly seen in Fig. 14 where the ebb current has its maximum near the land to the west whereas the flood current maximum is deflected farther to the east. The fact that ebb and flow current show different regional preferences makes it impossible to evaluate the Wadden Sea physics by only analyzing point measurements. Figure 14, however, shows the significant dynamic ebb–flood differences along the considered section. During flood, stratification is mostly unstable and eddy diffusivity is thus enhanced with the result that the fast-settling SPM is strongly mixed. During ebb the opposite happens: stratification is stable, eddy diffusivity is reduced, and SPM is much less mixed. Furthermore, during flood the current speed in the upper half of the water column is much more homogeneous during flood than during ebb. These dynamics compare well to the one-dimensional model results shown in Figs. 6 and 7.
Figure 15 shows the accumulation of SPM in the budget boxes I and II of the inner Sylt–Rømø Bight, which is done by comparing integrated SPM with total water volume in these budget boxes. When SPM is initialized to unity after 45 tidal periods, both values are identical. During the course of 10 further tidal periods, the total SPM for the fast-settling fraction is slowly increasing in comparison to the water volume, indicating an SPM accumulation in the inner Sylt–Rømø Bight. This is only the case for the simulation considering density gradients. The accumulation is smaller for the northern budget box II, the volume of which is much smaller than of the southern box I. If density gradients are neglected, an accumulation of the fast-settling SPM cannot be observed. For the slowly settling SPM, an accumulation in the inner Sylt–Rømø Bight is not visible for the simulation considering density gradients. For the simulation without density gradients, SPM is even slightly washed out of the inner bight, which can be explained by the net volume outflow due to the freshwater surplus.
Figure 16 finally shows the spatial pattern of the SPM accumulation. If density gradients are considered, SPM is significantly accumulating in the northeastern corner of the bight and in the southwestern part. Clear patterns of accumulation and depletion are also visible in the waters outside the bight, probably due to nonlinear interaction with the boundary conditions. When density gradients are not dynamically considered, hardly any SPM accumulation is visible, but some depletion is generated in the vicinity of the river runoff and on the sands due to net precipitation.
In this paper, we have shown that in general the water in the inner Wadden Sea is less dense than the outer water, with typical tidal density amplitudes of 0.5–1.0 kg m−3 (see Fig. 4). During winter, these density differences are mainly due to salinity differences, and during summer they are mainly due to temperature differences (see Fig. 4). This is due to differential effects of buoyancy fluxes (river runoff, net precipitation, surface warming), resulting in vertically averaged densities that are smaller for shallow water than for deep water.
Based on this characteristic feature, we have developed the hypothesis that density differences in the Wadden Sea contribute dynamically to the net SPM accumulation in the Wadden Sea. As mechanisms for this, we suggest the same ones as for the estuarine turbidity maxima that are typical for tidal estuaries (see section 2c). These mechanisms are based on small differences between ebb and flood that are difficult to quantify by means of field observations. However, because of the overlay of the permanent tidal activity with the generally existing along-current density gradients these small differences may result in long-term net onshore SPM fluxes.
Since mass-conserving numerical models with good process resolution have the capability to reproduce these long-term net fluxes, we rely on model simulations to challenge our hypothesis. The idealized three-dimensional simulations for the Sylt–Rømø Bight do indeed show a significant SPM accumulation in the inner bight if density gradients are considered and the settling velocity is fast enough. Without dynamic consideration of density differences no SPM accumulation is seen. This means that we did not find any evidence of net transport effects of settling lag, Stokes drift and other barotropic mechanisms relying on tidal asymmetries (see section 1), although we resolve these processes with the numerical model. Other possible processes such as scour lag and flocculation could not be quantified with our model since they are not included in the idealized SPM dynamics resolved here.
If the settling velocity of SPM is too small, accumulation is not seen in the model results also when density gradients are dynamically considered. The reason for this is that SPM with a settling velocity of only wc = 10−4 m s−1 sinks by only about 2 m during half a tidal cycle such that it is almost well mixed (see Fig. 8).
The three-dimensional model results are in good agreement with the one-dimensional results. The one-dimensional model predicts for the fast-settling SPM per tidal cycle a net onshore transport of 2.2% of the SPM transport amplitude. The three-dimensional model predicts about 10% net transport during 10 tidal cycles (see Fig. 15), which is less than for the one-dimensional model but of the same order of magnitude. The differences may be explained by the on average smaller water depth and the net outflow in the three-dimensional model. For slowly settling SPM, the one-dimensional model predicts a net transport of 0.33%, whereas the three-dimensional model with dynamic density differences predicts a zero net transport (see Fig. 15). Given that the net volume outflow out of the inner Sylt–Rømø Bight per tidal period is 0.1% of the mean volume and that again the mean water depth is less for the three-dimensional model, the SPM retention may result from a balance of a small net onshore transport and the net volume outflow. This supported by the decrease of SPM content in the inner bight by 0.2% per tidal period if density differences are not considered (see Fig. 15).
It is not only that we see a net onshore SPM transport into the inner bight for fast-settling SPM when dynamic effects of density gradients are considered, but that we also see the mechanisms directly acting when studying the dynamics along a cross section leading into the inner bight (see Fig. 14). During flood, stratification is unstable, thus increasing turbulent mixing with the consequence that velocity profiles are more homogeneous in the upper half of the water column and that SPM profiles are more mixed. During ebb, the opposite is seen.
The clearest signal of SPM accumulation is seen in the northern part of the inner bight (see Fig. 16). This is the location in the inner bight closest to an estuary: it is fairly narrow and shallow, with substantial freshwater supply, such that significant horizontal salinity gradients are resulting (see Fig. 11).
Of course, because of the idealized character of the simulations carried out here (no wind, no waves, periodic conditions, no SPM flux into sediment, no density differences due to SPM, no flocculation), and because of fairly coarse horizontal resolution (200 m), detailed SPM dynamics on small spatial scales cannot be reproduced. However, this numerical study strongly supports our hypothesis that dynamic effects of density gradients in the Wadden Sea significantly contribute to SPM accumulation. This is actually not in contrast to the general assumption that the Wadden Sea is vertically well mixed, since typical surface to bottom salinity differences are less than 0.1 psu (see Fig. 10).
From existing long-term observations these small vertical salinity gradients and the resulting small differences in ebb and flood salinity profiles cannot be significantly detected for various reasons. The piles installed in the Wadden Sea (see section 2a) are not installed in the deepest parts of the tidal channels in order to be outside of the shipping routes. The maximum vertical distance between salinity observations is 7.5 m for the Spiekeroog pile, which according to the numerical modeling results would amount to salinity differences of a few 1/100 psu. This is below the absolute accuracy of the CTD sensors used here (about 0.1 psu for the high-resolution measurements at the Spiekeroog pile) such that the discussed process cannot be assessed with this instrumentation. Even more substantial problems are present for the observation of current velocity profiles, where differences of about 3% between ebb and flood are typically overlaid by barotropic ebb to flood differences as shown in Fig. 12.
To further prove our hypothesis, targeted process observations in the field and detailed realistic numerical modeling studies will thus be necessary. Such observations need to include high-resolution SPM, current, density, and turbulence profiling along cross-sectional areas across tidal channels at various tidal phases.
With these methods it will be possible to estimate the dynamic impact of density gradients in a number of embayments in humid and tidally dominated regions such as the Jade Bay (German Bight), the Eastern Scheldt (Netherlands) and the Venice Lagoon (northern Adriatic).
The authors greatly appreciated the support of Karsten Bolding (Asperup, Denmark) for the maintenance of the public domain models GOTM and GETM. Furthermore, we are grateful for the critical comments on this manuscript by two anonymous reviewers.
Corresponding author address: Hans Burchard, Baltic Sea Research Institute Warnemünde, Seestraße 15, D-18119 Rostock, Germany. Email: firstname.lastname@example.org