The wave-induced solute flux from submerged sediment

The issue of the transport of dissolved nutrients and contaminants between the sediment in the bottom of a lake or reservoir and the body of water above it is an important one for many reasons. In particular the biological and chemical condition of the body of water is intricately linked to these mass transport processes. As the review by Boudreau (Rev Geophys 38(3):389–416, 2000) clearly demonstrates those transport processes are very complex involving mechanisms as diverse as the wave-induced flux between the sediment and the overlying water and the effect of burrowing animals on the transport within the sediment as well as basic diffusion mechanisms. The present paper focuses on one facet of these transport processes; we re-examine the balance of diffusion and wave-induced advection and demonstrate that the wave-induced flux of a solute from submerged sediment is not necessarily purely diffusive as suggested by Harrison et al. (J Geophys Res 88:7617–7622, 1983) but can be dominated by a mean or time-averaged flux induced by the advective fluid motion into and out of the sediment caused by the fluctuating pressure waves associated with wave motion. Indeed along the subtidal shoreline where the fluctuating bottom pressures are greatest, wave-induced advection will dominate the mean, time-averaged transport of solute into or out of the sediment as suggested in the work of Riedl et al. (Mar Biol 13:210–221, 1972). However, the present calculations also indicate that this advective flux decreases rapidly with increasing depth so that further away from the shoreline the advective flux becomes negligible relative to the diffusive flux and therefore the latter dominates in deeper water.


Introduction
This paper presents a re-evaluation of the wave-induced solute transport between lake-or ocean-bottom sediment and the overlying body of water. The issue is of particular importance to the biological and chemical condition of the overlying water because of the reactions that occur within the sediment and the associated transport of solute to and from the water. For example, nutrients generated by the decomposition of deposited organic mate-rial are recycled by transportation back up into the water. Also contaminants within the sediment may be similarly mixed up into the lake or ocean. For these reasons, significant attention has been given in the literature to these transport processes. In this paper we focus on the transport to and from the sediment that occurs in the subtidal zone. The transport that occurs along the beaches and in the tidal zone is much more complex involving factors such as breaking wave transport, tidal effects and the influence of freshwater runoff (see, for example, [1]). In subtidal waters, the possible transport mechanisms are fewer though still diverse: in addition to porous-media molecular diffusion they include wave-induced advection, tideinduced advection, the transport enhancement due to wave-induced ripples on the sediment surface, transport due to animal burrowing (see [2,3]) and that due to gas bubbles (or hot water) percolating up through the sediment [4]. Of these, this paper will be concerned only with wave-induced advection. Enhancement of transport by tide-induced advection has been detected in some North Sea sediments by [5] (purely diffusive transport was measured in other areas). The ripple-induced enhancement that has been measured in the laboratory by [6] and analysed by [7,8]. It probably pertains in some locations and is clearly related to the detailed investigations carried out by [9][10][11] and [12,13] on the solute flux between a stream and migrating sediment bedforms. Also site specific are the transport due to animal burrowing and percolating bubbles.
The focus of the present paper is restricted to the possible enhancement of the solute transport in the subtidal sediment caused by wave-induced pressures and flows. In particular, in the absence of the other possible transport mechanisms listed above, is there significant wave-induced advective enhancement of the solute transport or is porous-media diffusion the dominant process? This question has received significant attention in the literature, for example in the reviews by [2] and [14]. Transport in excess of that by molecular diffusion has certainly been observed by a number of researchers. Some have ascribed that excess to an increase in the effective diffusivity in the sediment due to tortuosity [15] or to turbulent mixing [16]. However experimental measurements of solute transport such as those of [17] off the coast of Georgia and Florida suggest that advective enhancement is important in some locations.
The literature contains a surprising diversity of analyses reflecting a significant divergence in what is considered the primary resistance to solute transport. Key analyses have been published by [18] and by [15] who present two different understandings and models of the wave-induced transport processes. The model of [15] contains no purely advective time-averaged transport due to flows into and out of the sediment that might be induced by the pressure fluctuations caused by surface or internal waves. Indeed the only effect that is included in their model is purely diffusive transport and the sole effect of the sediment is the modification of the porous medium diffusivity caused by the tortuosity of the diffusion pathway within the medium and by possible elevated levels of the diffusivity at higher Peclet numbers. In an alternative approach, [18] imply that the primary wave effect is due to advection into and out of the sediment. They associate the exchange of nutrients or contaminants between the sediment and the overlying water with the total volume of liquid that passes into and out of the sediment during the passage of one wavelength. Their analysis focuses on the total volume "filtered through the sediment" by adding the volume flowing into and out of the sediment (per unit surface area per unit time averaged over one cycle of the wave motion). [18] focus on this total advective transport volume; no diffusive contribution is included in their analysis. In this paper we will attempt to resolve this conflict and to identify those circumstances in which either of the two mechanisms, diffusion or advection, might dominate the time-averaged solute transport into or out of the sediment.

Transport equation
The basic equations for solute transport in a granular bed have been carefully developed and presented by [18] and [15] and need only be quoted here. The transport equation for the concentration, c, of a solute (which could be a contaminant or a nutrient dissolved in the liquid) is: where t is time, u is the pore water velocity field, and D is the diffusivity tensor. [15] address the circumstances in which the solute is being absorbed in the sediment, but only a minor modification is required to handle cases in which solute is produced rather than absorbed. To do so the term −qc in their version of the transport equation (where qc is the mass rate of absorption of solute per unit volume of the sediment per unit time) has been replaced in Eq. 1 by a mass rate of production per unit sediment volume per unit time of +q (c ∞ − c) where +q is the solute production factor and c ∞ is the concentration deep in the bed far from the sediment surface. We note that this form that includes the factor (c ∞ − c) is necessary in order for a steady state solute concentration to exist.
We also note as [15] observe that even in homogeneous isotropic sediments, the diffusivity in Eq. 1 must be represented by a tensor because of anisotropy caused by the direction of the pore water velocity vector, u. Experiments show (e.g. [19]) and analyses confirm [20] that in isotropic sediments the dispersion tensor in a coordinate system oriented with its first axis parallel to the pore water velocity vector u, has the following diagonal elements (the offdiagonal elements are zero): along this axis the "longitudinal" dispersion coefficient, denoted by D L , may be much larger than the "transverse" dispersion coefficient, D T , perpendicular to u. Both dispersion coefficients approach the porous media molecular diffusivity, D P , when the Peclet number, Pe = |u|d/D M 1 (d is some measure of the grain size, most conveniently the equivalent spherical diameter). Note also that D P is somewhat less than the pure liquid molecular diffusivity, D M , due to the more tortuous path for diffusion within the porous medium (usually D P ≈ 0.67D M ). More notably for Pe 1, D M > D T , usually by an order of magnitude. It is the combined effect on the diffusivity of the tortuosity and of the high Peclet number that [15] explore.
We define a coordinate system, (x, y, z), in which x and z are horizontal coordinates, y is vertical down and the pore velocity in the z direction is zero. Defining θ to be the angle of the pore velocity, u, relative to a horizontal plane [15] identify the following form of the diffusion term in the transport equation: where the subscripts x, y refer to partial differentiation with respect to those variables and α is the mean dispersion coefficient given by Note that when Pe 1, α → D P and β → 0.

Pore water velocity field
The pore water velocity may be evaluated using the incompressibility condition and Darcy's law for incompressible flow in a porous medium [21]. We focus here on the flows and solute fluxes induced by waves in the overlying liquid; both traveling surface waves and internal waves will be considered and it is convenient to consider an individual wavenumber, k, in the x direction associated with a radian frequency, ω, so that the pore velocity field, (u, v, w), can be represented by where t is time and U is the velocity fluctuation amplitude at the sediment surface. Note that the angle, θ , which the velocity vector makes with the x axis is given by θ = (kx − ωt).
According to the simplest form of Darcy's law, the sediment surface pressure variation, p(x, t), associated with the above velocity field is where ρ is the liquid density, g is the acceleration due to gravity and κ is the hydraulic conductivity of the sediment divided by the porosity. The velocity amplitude U is thus given by where P is the amplitude of the pressure fluctuations acting on the sediment surface as a result of the waves in the overlying water. The above applies to all traveling wave motions of small amplitude, whether surface waves or internal waves. In the particular case of linear surface waves of amplitude a on a homogeneous ocean of depth, H : but other forms of P and ω are readily incorporated into the system of governing equations including the effect of sediment surface slope [18,21]. The validity of such expressions for the fluid velocity within the sediment have been well tested by the extensive experimental measurements of [18] and others. [15] make note of additional limitations in the simple form of Darcy's law used in the above system of equations, including the neglect of fluid inertial terms (see also [22,23]). Finally, we note that the advective term in Eq. 1 can be written in terms of the above velocity field as u.∇c = U e −ky (c x cos θ − c y sin θ) (8) and the corresponding form of the solute transport Eq. 1 becomes To investigate the solute flux due to the above wave motion, we seek a solution of this equation which is purely periodic in x, t and therefore θ . Consequently, and, if we assume that α and β are constant and uniform and that there is no change in the z direction, then Eq. 9 becomes Though this appears complicated it has the advantage of being linear and homogeneous in the concentration, c. It is useful to rewrite it using a non-dimensional y coordinate, s, defined by Then, using the dimensionless parameters the governing Eq. 11 becomes is the governing equation to which we seek solutions. It can be greatly simplified by considering the typical magnitude of the various parameters as we will do in the sections which follow and in the Appendix. Once a solution is obtained we will evaluate the solute flux at the sediment surface in the downward direction (into the sediment), F(θ ), from and we will be most interested in the flux,F, averaged over one period of θ .

In the absence of waves
First, however, we note that in the absence of waves (k → 0 and ω → 0) the governing equation reduces to the simple diffusion equation to which the steady state solution is where c y=0 is the concentration at the sediment surface which would be uniform and constant in the absence of waves. It follows that the thickness of the "diffusion-reaction sublayer" just underneath the sediment surface (in which most of the solute production and diffusion occurs) is given by In this elementary case in which v y=0 = 0, the flux,F, simply becomes Note that the absence of a contribution to F(θ ) (Eq. 15) from the advective term is due to the imposed constant and uniform value of c y=0 . Equation 15 is also the result that [15] obtain when waves are present because they impose the condition that c is constant on the sediment surface. As a consequence, they are left only with the effects on the diffusivity of tortuosity and of high Peclet number. Thus, in the absence of waves and in the Harrison et al. wave effect analysis, the flux,F, divided by the flux,F M through a pure, non-moving liquid "sediment" with diffusivity D M , becomes

Simplification
We proceed with the evaluation of the effect of waves by examining some simplifications to the governing Eq. 14. The first step in this process is to gauge the magnitude of the parameter, . From the definition 13 and Eq. 18, is readily seen to be the ratio of the typical diffusionreaction sublayer thickness, (α/q ) 1 2 , to the wavelength divided by 2π, 1/k. The diffusionreaction sublayer thickness is usually of the order of millimeters to several centimeters [24,25] compared with wavelengths of the order of meters. Though much larger diffusion-reaction sublayer thicknesses are possible we will focus here on thicknesses that are smaller than meters. The parameter is therefore much less unity in many cases of practical interest and we will confine our attention to this range. Since η is always less than one, it follows that when 1 we need only be concerned with the sediment layer in which s < −1 , and the second and third lines of Eq. 14 become negligible (τ must also be small and this is addressed in Appendix). Equation 14 then takes the manageable, large wavelength form A detailed solution to this governing equation is delineated in Appendix. In that appendix we demonstrate that within some reasonable parametric limitations, an appropriate, approximate solution to the governing Eq. 21 which is consistent with the desired form of the sediment surface boundary condition (Eq. 24) leads to a time-averaged flux into the sediment,F, given byF wherec 0 andc are the time-averaged solute concentration and the fluctuating solute concentration at the sediment surface. Not surprisingly therefore, the boundary condition at the sediment surface plays a crucial role in determining the magnitude of the solute flux and we turn to that issue in the sections which follow. (We note that Eq. 22 could, of course, have been constructed heuristically without any of the complicated algebra presented in Appendix. However, the authors believe it is still useful to delineate the parametric limitations and assumptions behind that result.) The above expression, Eq. 22, displays an additional advective contribution, −c U/2, over and above that in the absence of waves or in the analysis of [15]. Consequently a key issue is the relative importance of diffusive and advective contributions to the time-averaged solute flux. From Eq. 22 the ratio of the advective flux to the diffusive flux is given by Therefore the issue depends on the product of two key quantities, the parameter τ and the ratio of the amplitude of the sediment surface concentration fluctuation to the difference between the deep sediment concentration and the time-averaged sediment surface concentration. We examine each of these quantities in the sections which follow, beginning with an estimate of τ followed by an examination of the sediment surface boundary condition and the quantitỹ c/2(c ∞ −c 0 ).

Parameter estimation
Measurements and estimates by [18,25] and others allow estimates of some of the controlling parameters, particularly the key parameter, τ = U/(q α) 1 2 . Beginning with the velocity U in the sediment, we will make use of the extensive measurements by [18] of the volume of liquid Q (per unit surface area per unit time averaged over a cycle of the wave motion) that passes into and out of the sediment during the passage of a wave. Note that Q = U/2π 2 in the present notation. Reidl et al. made measurements and estimates of Q for the sediment of the Atlantic Continental Shelf and concluded that Q varied from values of the order of Q = 2×10 −6 m/s (U = 4×10 −5 m/s) within a few miles of the coast to values smaller by an order of magnitude or more further out on the Continental Shelf. This decrease is primarily due to the increasing depth and its effect on the wave-induced bottom pressures. We will call this data for U and the associated dimensionless parameters, "coastal" data having U values of the order of 4 × 10 −5 m/s and decreasing with increasing depth.
However, sample calculations readily demonstrate that much larger values of U will pertain along the shoreline and, in the discussion, we will refer to these as "shoreline" data. For example, using Eq. 6 a 1 m amplitude wave of wavelength 10 m in shallow water with P/ρg = 1 m will produce U values of 10 −4 for the 130 µm Oklahoma 90 sand of [9] (which has a κ = 0.00011 m/s) and U = 10 −3 m/s for Elliott's 470 µm Ottawa sand (which has a κ = 0.0014 m/s). These much larger values are consistent with the measurements by [26] of the pore water velocities in the breaker zone of a laboratory beach. They imply very different solute transport along the shoreline.
Next we note that a typical value for the molecular diffusivity D P is 10 −9 m 2 s −1 [15] for fine sand, corresponding to a value of α ≈ 4 × 10 −9 m 2 s −1 . Moreover (α/q ) 1 2 is the typical diffusion-reaction sublayer thickness in the sediment and this thickness is estimated and measured to be about 0.03 m for fine sand [18,25]. It follows that with coastal U values of 4 × 10 −5 m/s and less and with fine sand the pertinent values of τ are 300 and less. On the other hand the above shoreline values of U lead to τ values of about 10 3 for the Oklahoma 90 sand and 10 4 for the Ottawa sand. In summary, the key parameter τ = U/(q α) 1 2 appears to decrease from values of the order of 10 4 for coarse sand (less for fine sand) at the subtidal shoreline to values of the order of 300 for fine sand in deeper coastal waters.
As a postscript we also note that other parameters of interest include , the ratio of the typical diffusion-reaction sublayer thickness to the wavelength (divided by 2π), which we have earlier concluded is much less than unity in almost all cases of practical interest and the parameter η which is clearly of order unity or smaller. We also note that, by definition, (q ) −1 is a measure of the time required for the source type used to deplete the concentration in the sediment by diffusion and since (ω) −1 is typical of the period of the wave motions we might reasonably expect that in almost all circumstances of interest the former time is much longer than the latter so that ξ 1 and ξη 1.

Sediment surface boundary
Clearly the time-averaged solute transport may depend crucially on the boundary condition(s) at the sediment surface. If the conditions are such that the solute concentration remains constant at the sediment surface (as [15] assume) then, sincec ≈ 0, the transport is purely diffusive; there would be no advective transport and the only significant effect on the timeaveraged transport would be the modification of the porous medium diffusivity caused by the tortuosity of the diffusion pathway within the sediment as suggested by [15]. This is because whenc = 0 the fluid entering and leaving the bed has the same concentration at all times and the net advective flux averaged over one cycle of the wave motion would be zero. In their alternative approach, [18] associate the exchange of nutrients or contaminants between the sediment and the overlying water with the total volume of liquid Q (per unit surface area per cycle of the wave motion), that passes into and out of the sediment during the passage of a wave (in the present notation Q = U/2π 2 ); no mention is made of any diffusive contribution. Moreover, the relation between that volume flux and the actual solute transport is not clarified.
In taking a closer look at the sediment surface concentration, we shall not assume a priori that c y=0 is constant but rather that it may change depending on whether the instantaneous advective flow is into or out of the sediment. We shall heuristically assume that c y=0 is correlated with the fluid velocity normal to the surface, namely v y=0 (positive downward) which we previously denoted by v y=0 = −U sin θ (note that for 0 < θ < π the flow is up out of the sediment whereas for π < θ < 2π the flow is into the sediment) and that c y=0 =c 0 +c sin θ (24) (Note that whenc is positive the surface concentration is greater when the flow is emerging from the sediment and less when the flow is into the sediment). Consider first the half-cycle during which flow and solute are emerging from the sediment into the overlying water in a lake, reservoir or ocean. After emerging it will be mixed by the turbulence in the benthic boundary layer [25]. During that half-period a layer of thickness, U/ω would be ejected from the sediment. Therefore, the extent to which the surface concentration will change during that half-period will depend on the ratio of U/ω to the thickness of the potential diffusion-reaction sublayer in the sediment, (α/q ) 1 2 . In other words it will depend on the parameter, U/ω(α/q ) 1 2 . If that parameter is much larger than unity then the change in the surface concentration will be comparable with the overall concentration difference, (c ∞ −c 0 ) since the diffusion-reaction sublayer in the sediment would be rapidly ejected and the surface concentration would approach c ∞ . On the other hand, if the parameter is much smaller than unity the surface concentration fluctuation will be much less.
Using the values for U and (α/q ) 1 2 from the preceding section and radian frequencies for the waves of ω = 0.6 → 6.0 s −1 leads to maximum coastal values of U/ω(α/q ) 1 2 of 2 × 10 −4 → 2 × 10 −3 and to shoreline values of U/ω(α/q ) 1 2 up to 0.5 for the coarse sand. We conclude that during the emerging flow half-cycle the surface concentration in the coastal zone and beyond will change very little but that the surface concentration in the shoreline zone will change significantly. Now consider the other half-period when the flow is into the sediment. Clearly a key consideration will be the mixing in the benthic boundary layer next to the bottom of oceans, lakes and reservoirs for that will determine the solute concentration entering the sediment. This mixing has been the subject of much research and debate (see, for example, [14,[27][28][29][30]). It can be crudely modelled using a turbulent diffusivity, D t , in the waters above the sediment. This allows an estimate of the potential diffusion layer thickness in the water in the absence of advection of (D t /ω) 1 2 (averaged over an integer number of wave periods). We seek the ratio of this to the thickness of layer ingested to the sediment during the same period, namely U/ω. If that ratio, (D t ω) 1 2 /U , is small then the diffusive benthic layer would be rapidly ingested and the surface concentration would approach the concentration in the water above that layer. If, on the other hand, that ratio is large the surface concentration fluctuation will be much less. Typical values of D t are the 3 × 10 −6 m 2 s −1 measured by [31] in Lake Kinnaret in Israel (see also [29]) and the 10 −5 m 2 s −1 measured by [25] in the Baltic Sea. Assuming D t ≈ 10 −5 m 2 s −1 it follows that for the coastal zone with U values of 4 × 10 −5 m s −1 (and smaller) and a radian wave frequency of ω = 0.6 → 6.0s −1 , the ratio, (D t ω) 1 2 /U , is in the range 60 → 200 though with much larger values for smaller U . On the other hand with subtidal shoreline values of U and with D t ≈ 10 −5 m 2 s −1 the parameter (D t ω) 1 2 /U could be as small as 0.2 → 8; however since the effective "turbulent" diffusivity in the shoreline zone is likely to be much greater than 10 −5 m 2 s −1 we estimate that the effective value for (D t ω) 1 2 /U at the shoreline is also likely to be very much greater than unity. These values suggests a mimimal change in the surface concentration during the sediment inflow half-period in both the coastal and shoreline zones.
We conclude from the preceding two paragraphs that the sediment surface concentration varies by very little in the coastal zone and beyond but may change significantly in the shoreline zone. Therefore, factoring in the respective τ values, we conclude that the ratio of the advective flux to the diffusive flux, given by Eq. 23 as τc/2(c ∞ −c 0 ), will be large in the shoreline zone (perhaps as much as 10 3 → 10 4 and therefore that solute transport at the shoreline will be dominated by wave-induced advection. On the other hand, the wave-induced advective transport in the coastal zone may be of the same order as the diffusive flux closest to the coast but will rapidly decrease with water depth so that, in deeper waters, the diffusive transport will dominate.
Parenthetically we note that, for simplicity, we have tacitly assumed similar grain sizes in the coastal and shoreline zones. When these differ substantially, the different grain sizes need to be included in evaluating the zonal sediment hydraulic conductivities and diffusivities and the relevant zonal parameters in which they appear.

Concluding remarks
We have demonstrated that the wave-induced flux of a solute from submerged sediment is not necessarily purely diffusive as suggested by [15] but can be dominated by a mean or time-averaged flux induced by the advective fluid motion into and out of the sediment caused by the fluctuating pressures associated with wave motion. For convenience we include the flowchart, Fig. 1, which illustrates the dominant processes and steps involved in evaluating and comparing the diffusive and advective contributions to the flux out of the sediment. Along the shoreline where the fluctuating bottom pressures are greatest, wave-induced advection will dominate the mean, time-averaged transport of solute into or out of the sediment. In this zone, the time-averaged transport will therefore be dominated by the mechanism suggested in the work of [18]. However, the present calculations also indicate that this advective flux decreases rapidly with increasing water depth so that further away from the shoreline the advective flux becomes negligible relative to the diffusive flux and, consequently, the diffusive flux dominates in deeper water where the calculations of [15] are therefore relevant.
Whether one can therefore assume that most of the solute transport occurs in the immediate vicinity of the subtidal shoreline or whether it includes significant contributions from deeper waters depends on the relative magnitudes of the concentration differences in the two zones. Clearly, models of the mixing processes in the entire body of the lake, reservoir or ocean and of the production/adsorption processes in the sediments at various depths will be needed to determine the prevailing concentration differences that drive mass transport. It is unlikely that these concentration differences would be similar in the subtidal shoreline and coastal zones. The present paper indicates the appropriate mechanism for transport in the different zones.