- Email: [email protected]
Modelling of Water Flows and Cohesive Sediment Fluxes in the Humber Estuary, UK
PII:
Marine Pollution Bulletin Vol. 37, Nos. 3±7, pp. 182±189, 1998 Ó 1999 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0025-326X/99 $ ± see front matter S0025-326X(99)00103-4
Modelling of Water Flows and Cohesive Sediment Fluxes in the Humber Estuary, UK Y. WU *, R. A. FALCONER and R. J. UNCLESà Cardi School of Engineering, Cardi University, P. O. Box 686, Cardi, Wales CF2 3TB, UK àCentre for Coastal Marine Science, Plymouth Marine Laboratory, Prospect Place, The Hoe, Plymouth PL1 3DH, UK
Details are given of a re®ned, three-dimensional, layerintegrated mathematical model for predicting accurately the ¯ow and cohesive sediment transport processes in estuarine and coastal waters. The model is applied to the Humber Estuary, UK, in order to predict the suspended cohesive sediment transport ¯uxes. Field-measured data at a station in the Trent Falls region (53° 42.80 N, 0° 37.80 W) and Station SG13 in the plume region (53° 35.190 N, 0° 13.900 E) were used to drive the model. Comparisons between the model predictions and ®eld-measured data at Station SG23 (53° 35.710 N, 0° 02.170 E) within the estuary (Hawke Anchorage) were made. Comparisons of the suspended cohesive sediment concentrations for the numerical model results and the ®eld data at Hawke Anchorage showed good agreement in the regime of deposition. In the regime of erosion, the model also gave reasonably good predictions, although the large variation in suspended sediment concentrations with depth was not well reproduced by the numerical model. Ó 1999 Elsevier Science Ltd. All rights reserved Keywords: estuarine and coastal hydrodynamics; cohesive sediment transport; deposition and erosion; numerical models; layer-integration; ULTIMATE QUICKEST; mass conservation.
Introduction The movement of cohesive sediments can cause siltation to occur in waterways and harbour docks, or erosion of estuarine banks. To deal with these problems, expensive dredging operations or bank protection are generally needed. The heavily contaminated materials resulting from industrial and municipal euents, or accidental oil spills, can also release heavy metals, mineral oils and other toxic contaminants, which may then be adsorbed onto ®ne sediments and be available for resuspension by strong tidal currents, short wave *Corresponding author.
182
action and dredging operations. This desorption of contaminants from their particulate phase can impact signi®cantly on the ecological balance of shallow, nearshore areas. Therefore, accurate predictions of cohesive sediment transport processes in estuarine and coastal waters is of vital importance for their environmental management. With this growing awareness of the problems associated with cohesive sediments, there has been considerable eort undertaken in recent years to develop and apply numerical models to predict coastal and estuarine cohesive sediment transport. For example, the two-layer model of Odd and Owen (1972); one-dimensional models, such as those by Scarlatos (1981) and Lin et al. (1983); two-dimensional, ®nite-element models based on solutions of the depth-averaged advective±diusion equation (e.g. Ariathurai and Krone, 1976; Hayter and Mehta, 1986); two-dimensional, ®nite-dierence models, based on solutions of the depth-averaged, advective± diusion equation (e.g. Cole and Miles, 1983; Falconer and Chen, 1996); and three-dimensional models (e.g. Markofsky et al., 1986; Nicholson and OÕconnor, 1986; Yost and Katopodes, 1996). In this paper a three-dimensional, layer-integrated numerical model has been re®ned to predict cohesive sediment transport ¯uxes in the Humber Estuary, UK. In the three-dimensional model, a split algorithm has been used to separate the three-dimensional transport equation into a horizontal, two-dimensional equation and a vertical, one-dimensional equation, because of the difference in length scales between the horizontal and vertical planes. A conservative form of the horizontal, twodimensional, layer-integrated advective±diusion equation has been used to ensure conservation of the cohesive sediment mass ¯uxes. A two-dimensional QUICKEST scheme has been developed by the authors, based on LeonardÕs (1979) one-dimensional QUICKEST scheme, with the scheme being used to solve the horizontal two-dimensional layer-integrated transport equation. The authors have also modi®ed the one-dimensional ULTIMATE algorithm, which has been used to avoid
Volume 37/Numbers 3±7/March±July 1998
unphysical numerical oscillations due to numerical dispersion (see Leonard, 1991; Wu and Falconer, 1998). For the vertical part of the advective±diusion equation, the Power-law implicit scheme has been used to predict accurately the solute distribution over depth. The model has been set up and applied to predict cohesive sediment transport ¯uxes in the Humber Estuary, from Trent Falls to the mouth of the estuary. Field-measured data at Trent Falls and Station SG13 were used to provide boundary conditions to drive the model, with comparisons between the model predictions and ®eld-measured data within the estuary also being undertaken.
The Mathematical Model Hydrodynamic model The hydrodynamic model used to predict the water elevations and velocity ®elds in coastal and estuarine waters initially involved solving the governing equations of ¯uid ¯ow, or the Navier±Stokes equations. Generally, the ¯ow was assumed to be isothermal and the vertical acceleration was assumed to be small compared to the gravitational acceleration, yielding a hydrostatic pressure distribution. Besides the gravity force in the vertical direction, the Coriolis force, as introduced by the earthÕs rotation, was also included in the horizontal plane, since the eects of the earthÕs rotation are often signi®cant in lateral circulation and mixing processes in estuarine and coastal ¯ows. The laminar shear stresses are generally small in comparison with the turbulent shear stresses and are therefore normally neglected. The turbulent shear stresses were equated to the wind stress on the free-surface and to the bottom friction on the bed, with both stresses being expressed using a quadratic friction law, including the local wind and bed layer velocities, respectively. Full details of the governing equations of ¯uid motion used in the model are given in Falconer (1993). A number of numerical models have been developed during the past two decades to predict hydrodynamic ¯ows in three-dimensions. For example, the layer-integrated models of Leendertse et al. (1973) and Lin and Falconer (1997), the ®nite-dierence, sigma coordinate transformation model of Shen (1987), and the ®nite-element model based on analytical, or semi-analytical, vertical velocity pro®les of Davies (1983). The layer-integrated method was used in this study, since models based on this method are more appropriate to treating the ¯ooding and drying processes, which are important in the Humber Estuary. Also, this modelling approach is computationally ecient for practical applications. In the layer-integrated model, the three-dimensional governing hydrodynamic equations were integrated within layers to give the layer-integrated equations. The layers and the relative variable locations in the x ÿ z plane are illustrated in Fig. 1. The horizontal viscosity was assumed to be constant in the vertical and equated
to the depth-averaged eddy viscosity in the current study (Falconer, 1993), with the vertical eddy viscosity ev being represented by a two-layer mixing length model (Rodi, 1984) of the following form: " 2 #1=2 2 ou ov ; 1 e v l2 oz oz where l is the mixing length, de®ned as: l jz
for jz 6 0:1H ;
l 0:1H
for jz > 0:1H
2
and j the von Karman's constant ( 0.41), and H the total water depth. The governing hydrodynamic equations were solved using a combined explicit and implicit ®nite-dierence scheme. Firstly, an alternating-direction implicit scheme, based on FalconerÕs (1986) model, was used to solve the depth-integrated hydrodynamic equations to give the water elevation ®eld. The layer-integrated equations were then solved to obtain the layer-integrated velocities, using the water elevation ®eld predicted by the depth-integrated equations. In solving the depth-integrated hydrodynamic equations, the layer-integrated velocities were integrated to obtain the depthmean velocity. The Crank±Nicolson scheme was used to solve the layer-integrated hydrodynamic equations, with the vertical diusion terms being treated implicitly and the remaining terms being treated explicitly. Two iterations were performed to solve the coupled, depth-integrated and layer-integrated equations. The ¯ooding and drying processes were modelled using a robust scheme developed by Falconer and Chen (1991). For details of the hydrodynamic model and the numerical scheme, see Lin and Falconer (1997). Sediment transport model The following advective±diusion, partial dierential equation was used to describe the suspended sediment transport processes: o/ o o o u/ v/ w ÿ ws / ot ox oy oz o o/ o o/ o o/ ex ÿ ey ÿ ez 0; ÿ ox ox oy oy oz oz
3
where / is the cohesive sediment concentration, u, v and w are ¯uid velocity components in the x, y and z directions, respectively, ws is the apparent sediment settling velocity, and ex , ey and ez are turbulent diusion coecients in the x, y and z directions, respectively. As a result of ¯oc aggregation due to inter-particle collisions and surface electro-chemical forces, cohesive sediments settle by ¯ocs rather than by individual particles. It was found that the settling velocity of ¯ocs depended strongly upon the suspended cohesive sediment concentration. The dependence of the settling velocity on the local concentration generally falls within one of the three following ranges (see Mehta, 1993): 183
Marine Pollution Bulletin
Fig. 1 Sketch of layers.
(i) Free settling / < /1 0:1 ÿ 0:3 g=l ws
1gD2s
s ÿ 18m
;
4
where m is the kinematic viscosity for clear water, s the speci®c gravity of suspended sediment and Ds the ¯oc diameter. (ii) Flocculation settling /1 < / < /2 0:3 ÿ 10 g=l ws k1 /4=3 ;
5
where k1 is an empirical coecient. (iii) Hindered settling / > /2 4:66
ws ws0 1 ÿ k2 / ÿ /2
;
6
where ws0 is the value of ws at the concentration level /2 , and k2 the inverse of /s ÿ /2 , where /s is the concentration at which ws 0. In estuarine and coastal waters, the horizontal length scale is generally much larger than the vertical scale, and hence an operator splitting algorithm may be used to separate the three-dimensional, advective±diusion equation into a horizontal, two-dimensional equation and a vertical, one-dimensional equation. Firstly, the two-dimensional advective±diusion equation was solved horizontally, and secondly, the one-dimensional vertical advective±diusion equation was solved vertically. To be consistent with the three-dimensional layerintegrated hydrodynamic model, the horizontal twodimensional equation was also integrated through the layers to give a layer-integrated two-dimensional equation. The two-dimensional QUICKEST scheme extended from LeonardÕs (1979) one-dimensional QUICKEST scheme was used to solve the horizontal two-dimensional layer-integrated equation. The modi®ed one-dimensional ULTIMATE algorithm was also used to prevent non-physical numerical oscillations, which frequently arise due to numerical dispersion when large concentration gradients exist (see Wu and Falconer, 184
1998). Another feature of the present model is that an additional source term associated with layer-integration of the free-surface ¯ow was introduced into the model, which was found to be of vital importance for mass conservation of pollutants (see Wu and Falconer, 1998). Boundary conditions The open boundary conditions were speci®ed by prescribing concentration data from ®eld measurements when in¯ow conditions prevailed, and with the concentration pro®les along the boundary being obtained by extrapolation using a ®rst-order upwind dierence scheme for out¯ow conditions. At bank (or closed) boundaries, normal derivatives of the concentration were set to zero. The one-dimensional, advective±diusion equation was solved using a non-uniform grid spacing in the vertical direction (see Fig. 1). Since the diusion term was very important for the vertical concentration distribution, and since some of the grid sizes were very small near the sea bed and water surface, the implicit power-law scheme of Patanker (1988) was adapted to avoid the use of very small time steps. For the sediment deposition and erosion rates, unlike two-dimensional depth-integrated models where these processes are modelled using sink and source terms, these ¯uxes are represented in the three-dimensional model using the following bed conditions: ÿ ws / ÿ ez
o/ qdep oz
when sb 6 sd
deposition; 7
ÿ ws / ÿ ez
o/ qero oz
when sb P se
erosion;
ÿ ws / ÿ ez
o/ 0 oz
when sd < sb < se
8
equilibrium; 9
where sb is the bed shear stress, sd the critical shear stress beyond which there is no further deposition, se the
Volume 37/Numbers 3±7/March±July 1998
critical shear stress for erosion, and qdep and qero represent the deposition and erosion rates, respectively, at the bed. In recent years, considerable eort has been made to study the mechanisms of cohesive sediment transport, and many experimental and ®eld studies have been carried out to investigate deposition and erosion rates of cohesive sediments (e.g. Krone, 1962; Mehta, 1973; Thorn and Parsons, 1977; Mehta, 1988; Parchure and Mehta, 1985; Thorn, 1981; Lick et al., 1995). The deposition rate proposed by Krone (1962) was used in this study: sb ; 10 qdep ÿws /b 1 ÿ sd where /b is the near-bed cohesive sediment concentration, with typical values of the critical shear stress sb for 2 deposition being 0.04±0.15 N=m (Krone, 1962; Mehta, 1973). Likewise, the erosion rate for soft natural mud can be represented by the following empirical expression (Raudkivi, 1990): qero 1=2 exp a sb ÿ se ; 11 qf where the empirical coecients are qf 0:0000042 and a 8:3 as suggested by Thorn and Parsons (1977). Typical values of the critical shear stress se for the ero2 sion of soft mud are 0.07±0.17 N=m . To complete the solution of the vertical, one-dimensional cohesive sediment transport equation, the following condition was applied at the free surface: w ÿ ws / ÿ ez
o/ 0 oz
12
which describes the condition that there can be no ¯ux across the free surface.
Humber Estuary Study General The Humber Estuary is one of the main estuaries in the UK, situated along the coast of northern England, and strategically important in terms of its ¯ow and solute input to the North Sea and its shipping links to mainland Europe. The main part of the estuary is about 62 km long, from Trent Falls to Spurn Head (see Fig. 2), and the tidal in¯uence extends for another 62 km up the River Ouse and about 72 km up the River Trent. As part of this study within the LOIS project, a comprehensive ®eld observation programme was undertaken and the model outlined herein was set up and applied to simulate cohesive sediment transport processes in the estuary from Trent Falls to just beyond Spurn Head. Field observations Observations of water elevations, currents, salinity, temperature and suspended particulate matter concentrations were made in the Humber mouth and plume
region of freshwater in¯uence during June 1995. Four tidal-cycle anchor stations were worked from R.V. Challenger, three of which were in the plume region: Stations SG13, SG24 and SG10. Station SG13 (53° 35.190 N, 0° 13.900 E) was occupied on 3 June between 05:17 and 18:45 h. The tidal range at Immingham was 4.9 m (mean tides). Station SG24 (53° 32.060 N, 0° 15.290 E) was worked on 4 June between 05:47 and 18:47 h. The tidal range was 4.5 m (large neaps). SG10 (53° 36.820 N, 0° 11.350 E) was worked on 7 June between 06:05 and 19:02 h. The tidal range was 3.8 m (neaps). The bed sediments at these stations were broken shells, coarse sand and larger material. The fourth Station SG23 (53° 35.710 N, 0° 2.170 E) was positioned within the estuary at Hawke Anchorage and was worked on 8 June between 07:21 and 20:53 h. The tidal range was 3.9 m (neaps). The bed sediments at this site were ®ne sand and smaller material. A tidal-cycle anchor station was worked from R.V. Sea Vigil, in the Trent Falls region of the upper Humber (53° 42.80 N, 0° 37.80 W). The station was occupied on 18 July between 09:04 and 21:54 h. The tidal range at Immingham was 5.3 m (small springs). The bed sediments were very ®ne sand and mud. The experimental procedure comprised water-column pro®ling from surface to bed at 1 or 2 m depth intervals, dependent upon the total water depth, and at approximately 0.5 h time intervals during each tidal cycle. The ship remained at anchor throughout the observations. A direct reading current meter was used to measure the current speed and direction. Direct reading temperature, conductivity, salinity and depth sensors were deployed, together with multi-range transmissometers for the determination of suspended particulate matter concentrations. Instruments were calibrated by comparison with data from laboratory analyses of water sampled during the measurements. Measurements in the plume area, just outside of the estuary and the model boundary, showed that suspended particulate matter, SPM, concentrations were small. At SG13, for example, the mean and standard deviation (SD) of the SPM for the whole tide was 13 4 mg lÿ1 . The mean and SD of the salinity was 33.8 0.7. There was a strong inverse relationship between SPM and salinity that indicated both that the dominant source of SPM was the estuary and that SPM transport was essentially conservative on the intra-tidal time-scale. Similar results were obtained at SG24 and SG10. The situation was dierent within the estuary at SG23. Here the SPM concentrations were signi®cantly greater and there was evidence of resuspension of sediments at the strongest current speeds. For the whole tide, SPM and salinity were 46 25 mg lÿ1 and 31.4 1.3, respectively. Measurements at Trent Falls showed that SPM concentrations were much greater in the upper reaches. For the whole tide, the means and SD of SPM and salinity were 2.4 2.6 g lÿ1 and 14.5 4.2, respectively. There was a signi®cant inverse correlation with salinity that 185
Marine Pollution Bulletin
Fig. 2 Map of the Humber Estuary and the locations of ®eld observation stations.
indicated that much greater SPM concentrations existed within the Ouse and Trent sub-estuaries. These turbidity maximum regions of extremely high SPM concentrations comprised suspended ®ne sediment (<63 lm), in which the majority of sediment existed as ¯ocs of clay and very ®ne silt particles. Model application and results The model area was represented horizontally using a mesh of 118 56 uniform grid squares, with a grid size of 500 m. Vertically, eight layers were used, with the thickness of the top layer being 4 m at mean water level, and with the other layers each being 3 m thick. The streamline boundary was located such that the boundary was as parallel as possible to the dividing streamline that separated the north±south alongshore current and the estuarine in¯ow±out¯ow. The water elevation recorded at Station SG13 was chosen as the seaward boundary condition to drive the tidal current. Due to the lack of volume ¯ux data at Trent Falls (Fig. 2), a water elevation boundary condition was prescribed at this landward boundary. A bed roughness height of 10 mm was assumed, as suggested by model calibration for a previous study (see Lin and Falconer, 1997). To achieve better predictions of the cohesive sediment distributions, accurate modelling of the hydrodynamic conditions in the estuary was extremely important. Since measurements taken at Trent Falls and Station SG13 were not simultaneous, i.e. the measurements at Station SG13 were taken for mean tide, whereas conditions at Trent Falls were taken during the spring tide, the water elevations at Trent Falls had to be adjusted to mean tide conditions. Calibration was carried out during the adjustment according to the water elevations and mean velocities measured at Hawke Anchorage. A comparison of the velocities for the dierent layers is given in Fig. 3. As can be seen from the results, the 186
general variation in time in the measured and predicted velocities at Hawke Anchorage shows a reasonable degree of agreement with regard to the depth-mean velocities. However, the variation over depth of the measured horizontal velocity components is generally much greater than the variation predicted from the numerical model. This under-prediction of the variation in the horizontal velocity ®eld is primarily thought to be due to the inclusion of a relatively simple turbulence model (i.e. the two-layer mixing length model) and the relative coarseness of the vertical grid, particularly near the bed boundary. Furthermore, the vertical advectionmomentum terms are represented using only a quasisecond order accurate scheme. These simpli®cations are currently being investigated further through another research programme. Suspended sediment concentrations are relatively low at Hawke Anchorage so that the free settling criterion is applicable [Eq. (4)]. Due to a lack of information on the primary particle sediment sizes at Hawke Anchorage, a sediment-speci®c gravity of 2.65 and a grain (micro¯oc) size of 20 lm was used in Eq. (4). Larger aggregates may be present, but these will have a lower speci®c gravity and settle relatively slowly. The critical shear stress for 2 deposition was set to: sd 0:07 N=m , and the critical 2 shear stress for erosion assumed to be: se 0:15 N=m . These values were in the ranges proposed respectively by Krone (1962) and Thorn and Parsons (1977), and were re®ned by trial and error. Fig. 4 shows a comparison of the ®eld-measured data and the model-predicted suspended cohesive sediment concentration variations with time. As the concentration variations with depth from the model predictions were small, only the depth mean concentration is shown in Fig. 4. The results generally show reasonable agreement between the model predictions and the ®eld-measured data in the regime of deposition. However, in the regime of erosion, the ®eld data showed a relatively large
Volume 37/Numbers 3±7/March±July 1998
Fig. 3 Comparison of predicted and ®eld measured velocities at Hawke anchorage, taken on 8 June 1995.
Fig. 4 Comparison of predicted and ®eld measured cohesive sediment concentrations at Hawke anchorage, taken on 8 June 1995.
variation with depth that has not been reproduced in the model predictions. This discrepancy is primarily thought to be due to the exclusion of non-cohesive sediment transport predictions in the model simulations reported herein. For the case of cohesive sediment transport predictions, once the erosive shear stress has been exceeded then the ®ne particulate matter will rapidly be distributed fairly uniformly through the water column in a strongly mixed estuary, because the settling velocity is very small. In contrast, for the case of non-cohesive sediment transport predictions, the settling velocity is much greater and the vertical variation much more pronounced. For example, previous simulations of the vertical variation of non-cohesive suspended sediment ¯uxes in the Humber Estuary showed very good agreement with ®eld data through the water column at Sunk Channel and Middle Shoal (Fig. 2; see Lin and Falconer, 1996). Hence, the disparity between the measured
and predicted sediment concentration distributions through the water column is thought to be due to the existence of a mixture of cohesive and non-cohesive sediment particles occurring in-situ during the erosion phase, whereas only cohesive sediment predictions have been included in the numerical model. This degree of variation between the predicted and measured suspended sediment concentration distributions is currently being investigated further in a similar model study for the Mersey Estuary, where the modelling of sediment mixtures is being considered.
Discussion and Conclusions A re®ned, three-dimensional, layer-integrated numerical model has been applied to simulate cohesive sediment transport processes in the Humber Estuary, UK. A combined layer-integrated and depth-integrated 187
Marine Pollution Bulletin
scheme was used to solve the hydrodynamic equations. The original three-dimensional advective±diusion equation was split into a two-dimensional horizontal and a one-dimensional vertical equation to implement the variation between the horizontal and vertical scales. A split algorithm was used to separate the three-dimensional transport equation into a horizontal, twodimensional equation and a vertical, one-dimensional equation due to the dierent length scales in the horizontal and vertical planes. To be consistent with the three-dimensional layer-integrated hydrodynamic model, the two-dimensional horizontal advective±diusion equation was integrated over the layers to give the layerintegrated advective±diusion equation. A conservative form of the horizontal, two-dimensional, layer-integrated advective±diusion equation was used to ensure mass conservation of the cohesive sediment ¯uxes. A two-dimensional modi®ed QUICKEST scheme was developed and used in this study to achieve high numerical accuracy. To avoid numerical oscillations, a modi®ed one-dimensional universal limiter ULTIMATE scheme was also used. Comparisons between the measured and predicted velocities at Hawke Anchorage generally showed good agreement between the depth-mean velocities. However, the measured velocities showed a greater variation about the mean value through the water column. This was thought to be primarily due to the relatively simple turbulence model and the coarse grid resolution near the bed. The model is currently being extended to include a version of the two-equation k±e turbulence model in the vertical. Comparisons of the suspended cohesive sediment concentrations were also undertaken at Hawke Anchorage and showed good agreement for the deposition phase. However, during the erosion phase the depthmean concentrations were generally in good agreement, but again the variation over the water column was much greater for the measured data than for the predicted results. This disparity in the vertical concentration distribution during the erosion phase was primarily thought to be due to mixtures of cohesive and non-cohesive sediment ¯uxes occurring in-situ, whereas only cohesive sediment transport predictions were included in the model. These results have indicated that further research needs to be undertaken on predicting ¯uxes of sediment mixtures. This is LOIS publication number 551 of the LOIS Thematic Research Programme, carried out under a Special Topic Award from the Natural Environment Research Council. The authors are grateful to NERC for supporting this study. Ariathurai, R. and Krone, R. B. (1976) Finite element model for cohesive sediment transport. Journal of Hydraulic Engineering Division, ASCE 102, 323±338. Cole, P. and Miles, G. V. (1983) Two-dimensional model of mud transport. Journal of Hydraulic Engineering Division, ASCE 109, 1± 12.
188
Davies, A. M. (1983) Formulation of a linear three-dimensional hydrodynamic sea model using a Galerkin-eigenfunction method. International Journal for Numerical Methods in Fluids 3, 33±60. Falconer, R. A. and Chen, Y. (1991) An improved representation of ¯ooding and drying and wind stress eects in a 2-d numerical model. Proceedings of the Institution of Civil Engineers, Part 2, Research and Theory 91, 659±687. Falconer, R. A. (1993) An introduction to nearly horizontal ¯ows. In Coastal, Estuarial and Harbour Engineers' Reference Book, ed. M. B. Abbott and W. A. Price, pp. 27±36. E & FN Spon Ltd., London. Falconer, R. A. and Chen, Y. (1996) Modelling sediment transport and water quality processes on tidal ¯oodplains. In Floodplain Processes, eds. M. G. Anderson, D. E. Walling and P. D. Bates, pp. 361±398. Wiley, New York. Hayter, E. J. and Mehta, A. J. (1986) Modelling cohesive sediment transport in estuarine waters. Applied Mathematical Modelling 10, 294±303. Krone, R. B. (1962) Flume studies of the transport of sediment in estuarial processes. Final Report, Hydraulic Engineering Laboratory and Sanitary Engineering Research Laboratory. University of California, Berkeley. Leendertse, J. J., Alexander, R. C. and Liu, S. K. (1973) A threedimensional model for estuaries and coastal seas: Volume I. Principles of Computation. Report No. R-1417-OWRR, Rand Corporation, Santa Monica. Lick, W., Xu, Y. J. and McNeill, J. (1995) Resuspension properties of sediments from the Fox, Sagina and Bualo rivers. Journal of Great Lakes Research 21, 257±274. Lin, B. and Falconer, R. A. (1996) Numerical modelling of threedimensional suspended sediment for estuarine and coastal waters. Journal of Hydraulic Research, IAHR 34, 435±456. Lin, B. and Falconer, R. A. (1997) Three-dimensional layer-integrated modelling of estuarine ¯ows with ¯ooding and drying. Estuarine, Coastal and Shelf Science 44, 737±751. Lin, P., Huan, J. and Li, X. (1983) Unsteady transport of Suspended load at small Concentrations. Journal of Hydraulic Engineering, ASCE 109, 86±98. Markofsky, M., Lang, G. and Schubert, R. (1986) Suspended sediment transport in rivers and estuaries. In Lecture Notes on Coastal Engineering Studies, Physics of Shallow Estuaries and Bays, ed. J. Van de Kreeke, pp. 210±227. Mehta, A. J. (1973) Deposition behaviour of cohesive sediment, Ph.D. Thesis, University of Florida, Gainesville. Mehta, A. J. (1988) Laboratory studies on cohesive sediment deposition and erosion. In Physical Processes in Estuaries, ed. J. Dronkers and W. van Leussen, pp. 427±445. Springer, New York. Mehta, A. J. (1993) Hydraulic behaviour of ®ne sediment. In Coastal, Estuarial and Harbour Engineers' Reference Book, ed. M. B. Abbott and W. A. Price, pp. 577±584. E & FN Spon Ltd., London. Nicholson, J. and OÕConnor, B. A. (1986) Cohesive sediment transport model. Journal of Hydraulic Engineering, ASCE 112, 621±640. Odd, N. V. M. and Owen, M. W. (1972) A two-layer model for mud transport in the Thames Estuary. Proceedings of the Institution of Civil Engineers, London, Supplement (ix), pp. 175±205. Parchure, T. M. and Mehta, A. J. (1985) Erosion of soft cohesive sediment deposits. Journal of Hydraulic Engineering, ASCE 111, 1308±1326. Raudkivi, A. J. (1990) Loose Boundary Hydraulics. 3rd edn. Pergamon Press, Oxford. Rodi, W. (1984) Turbulence Models and their Application in Hydraulics. 2nd edn. IAHR Publication, Delft. Scarlatos, P. D. (1981) On the numerical modelling of cohesive sediment transport. Journal of Hydraulic Research, IAHR 19, 61± 68. Shen, Y. P. (1987) On modelling three-dimensional estuarine and marine hydrodynamics. In Three-Dimensional Models of Marine and Estuarine Dynamics, ed. J. C. J. Nihoul and B. M. Jamart, pp. 35± 44. Elsevier Oceanography Series, Elsevier, Amsterdam. Thorn, M. F. C. and Parsons, J. G. (1977) Properties of Grangemouth Mud. Hydraulic Research Station, Rep. No. EX 781, Wallingford, UK. Thorn, M. F. C. (1981) Physical processes of siltation in channels. Proceedings of the Conference on Hydraulic Modelling Applied to Maritime Engineering Problems, pp. 47±55. Institution of Civil Engineers, London. Wu, Y. and Falconer, R. A. (1998) Two-dimensional ULTIMATE QUICKEST scheme for pollutant transport in free-surface ¯ows.
Volume 37/Numbers 3±7/March±July 1998 Third International Conference on Hydroscience and Engineering. Cottbus/Berlin, Germany.
Yost, S. A. and Katopodes, N. D. (1996) 3D numerical modeling of cohesive sediments. Proceedings of the International Conference on Estuarine and Coastal Modeling, pp. 478±489.
189
Marine Pollution Bulletin Vol. 37, Nos. 3±7, pp. 182±189, 1998 Ó 1999 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0025-326X/99 $ ± see front matter S0025-326X(99)00103-4
Modelling of Water Flows and Cohesive Sediment Fluxes in the Humber Estuary, UK Y. WU *, R. A. FALCONER and R. J. UNCLESà Cardi School of Engineering, Cardi University, P. O. Box 686, Cardi, Wales CF2 3TB, UK àCentre for Coastal Marine Science, Plymouth Marine Laboratory, Prospect Place, The Hoe, Plymouth PL1 3DH, UK
Details are given of a re®ned, three-dimensional, layerintegrated mathematical model for predicting accurately the ¯ow and cohesive sediment transport processes in estuarine and coastal waters. The model is applied to the Humber Estuary, UK, in order to predict the suspended cohesive sediment transport ¯uxes. Field-measured data at a station in the Trent Falls region (53° 42.80 N, 0° 37.80 W) and Station SG13 in the plume region (53° 35.190 N, 0° 13.900 E) were used to drive the model. Comparisons between the model predictions and ®eld-measured data at Station SG23 (53° 35.710 N, 0° 02.170 E) within the estuary (Hawke Anchorage) were made. Comparisons of the suspended cohesive sediment concentrations for the numerical model results and the ®eld data at Hawke Anchorage showed good agreement in the regime of deposition. In the regime of erosion, the model also gave reasonably good predictions, although the large variation in suspended sediment concentrations with depth was not well reproduced by the numerical model. Ó 1999 Elsevier Science Ltd. All rights reserved Keywords: estuarine and coastal hydrodynamics; cohesive sediment transport; deposition and erosion; numerical models; layer-integration; ULTIMATE QUICKEST; mass conservation.
Introduction The movement of cohesive sediments can cause siltation to occur in waterways and harbour docks, or erosion of estuarine banks. To deal with these problems, expensive dredging operations or bank protection are generally needed. The heavily contaminated materials resulting from industrial and municipal euents, or accidental oil spills, can also release heavy metals, mineral oils and other toxic contaminants, which may then be adsorbed onto ®ne sediments and be available for resuspension by strong tidal currents, short wave *Corresponding author.
182
action and dredging operations. This desorption of contaminants from their particulate phase can impact signi®cantly on the ecological balance of shallow, nearshore areas. Therefore, accurate predictions of cohesive sediment transport processes in estuarine and coastal waters is of vital importance for their environmental management. With this growing awareness of the problems associated with cohesive sediments, there has been considerable eort undertaken in recent years to develop and apply numerical models to predict coastal and estuarine cohesive sediment transport. For example, the two-layer model of Odd and Owen (1972); one-dimensional models, such as those by Scarlatos (1981) and Lin et al. (1983); two-dimensional, ®nite-element models based on solutions of the depth-averaged advective±diusion equation (e.g. Ariathurai and Krone, 1976; Hayter and Mehta, 1986); two-dimensional, ®nite-dierence models, based on solutions of the depth-averaged, advective± diusion equation (e.g. Cole and Miles, 1983; Falconer and Chen, 1996); and three-dimensional models (e.g. Markofsky et al., 1986; Nicholson and OÕconnor, 1986; Yost and Katopodes, 1996). In this paper a three-dimensional, layer-integrated numerical model has been re®ned to predict cohesive sediment transport ¯uxes in the Humber Estuary, UK. In the three-dimensional model, a split algorithm has been used to separate the three-dimensional transport equation into a horizontal, two-dimensional equation and a vertical, one-dimensional equation, because of the difference in length scales between the horizontal and vertical planes. A conservative form of the horizontal, twodimensional, layer-integrated advective±diusion equation has been used to ensure conservation of the cohesive sediment mass ¯uxes. A two-dimensional QUICKEST scheme has been developed by the authors, based on LeonardÕs (1979) one-dimensional QUICKEST scheme, with the scheme being used to solve the horizontal two-dimensional layer-integrated transport equation. The authors have also modi®ed the one-dimensional ULTIMATE algorithm, which has been used to avoid
Volume 37/Numbers 3±7/March±July 1998
unphysical numerical oscillations due to numerical dispersion (see Leonard, 1991; Wu and Falconer, 1998). For the vertical part of the advective±diusion equation, the Power-law implicit scheme has been used to predict accurately the solute distribution over depth. The model has been set up and applied to predict cohesive sediment transport ¯uxes in the Humber Estuary, from Trent Falls to the mouth of the estuary. Field-measured data at Trent Falls and Station SG13 were used to provide boundary conditions to drive the model, with comparisons between the model predictions and ®eld-measured data within the estuary also being undertaken.
The Mathematical Model Hydrodynamic model The hydrodynamic model used to predict the water elevations and velocity ®elds in coastal and estuarine waters initially involved solving the governing equations of ¯uid ¯ow, or the Navier±Stokes equations. Generally, the ¯ow was assumed to be isothermal and the vertical acceleration was assumed to be small compared to the gravitational acceleration, yielding a hydrostatic pressure distribution. Besides the gravity force in the vertical direction, the Coriolis force, as introduced by the earthÕs rotation, was also included in the horizontal plane, since the eects of the earthÕs rotation are often signi®cant in lateral circulation and mixing processes in estuarine and coastal ¯ows. The laminar shear stresses are generally small in comparison with the turbulent shear stresses and are therefore normally neglected. The turbulent shear stresses were equated to the wind stress on the free-surface and to the bottom friction on the bed, with both stresses being expressed using a quadratic friction law, including the local wind and bed layer velocities, respectively. Full details of the governing equations of ¯uid motion used in the model are given in Falconer (1993). A number of numerical models have been developed during the past two decades to predict hydrodynamic ¯ows in three-dimensions. For example, the layer-integrated models of Leendertse et al. (1973) and Lin and Falconer (1997), the ®nite-dierence, sigma coordinate transformation model of Shen (1987), and the ®nite-element model based on analytical, or semi-analytical, vertical velocity pro®les of Davies (1983). The layer-integrated method was used in this study, since models based on this method are more appropriate to treating the ¯ooding and drying processes, which are important in the Humber Estuary. Also, this modelling approach is computationally ecient for practical applications. In the layer-integrated model, the three-dimensional governing hydrodynamic equations were integrated within layers to give the layer-integrated equations. The layers and the relative variable locations in the x ÿ z plane are illustrated in Fig. 1. The horizontal viscosity was assumed to be constant in the vertical and equated
to the depth-averaged eddy viscosity in the current study (Falconer, 1993), with the vertical eddy viscosity ev being represented by a two-layer mixing length model (Rodi, 1984) of the following form: " 2 #1=2 2 ou ov ; 1 e v l2 oz oz where l is the mixing length, de®ned as: l jz
for jz 6 0:1H ;
l 0:1H
for jz > 0:1H
2
and j the von Karman's constant ( 0.41), and H the total water depth. The governing hydrodynamic equations were solved using a combined explicit and implicit ®nite-dierence scheme. Firstly, an alternating-direction implicit scheme, based on FalconerÕs (1986) model, was used to solve the depth-integrated hydrodynamic equations to give the water elevation ®eld. The layer-integrated equations were then solved to obtain the layer-integrated velocities, using the water elevation ®eld predicted by the depth-integrated equations. In solving the depth-integrated hydrodynamic equations, the layer-integrated velocities were integrated to obtain the depthmean velocity. The Crank±Nicolson scheme was used to solve the layer-integrated hydrodynamic equations, with the vertical diusion terms being treated implicitly and the remaining terms being treated explicitly. Two iterations were performed to solve the coupled, depth-integrated and layer-integrated equations. The ¯ooding and drying processes were modelled using a robust scheme developed by Falconer and Chen (1991). For details of the hydrodynamic model and the numerical scheme, see Lin and Falconer (1997). Sediment transport model The following advective±diusion, partial dierential equation was used to describe the suspended sediment transport processes: o/ o o o u/ v/ w ÿ ws / ot ox oy oz o o/ o o/ o o/ ex ÿ ey ÿ ez 0; ÿ ox ox oy oy oz oz
3
where / is the cohesive sediment concentration, u, v and w are ¯uid velocity components in the x, y and z directions, respectively, ws is the apparent sediment settling velocity, and ex , ey and ez are turbulent diusion coecients in the x, y and z directions, respectively. As a result of ¯oc aggregation due to inter-particle collisions and surface electro-chemical forces, cohesive sediments settle by ¯ocs rather than by individual particles. It was found that the settling velocity of ¯ocs depended strongly upon the suspended cohesive sediment concentration. The dependence of the settling velocity on the local concentration generally falls within one of the three following ranges (see Mehta, 1993): 183
Marine Pollution Bulletin
Fig. 1 Sketch of layers.
(i) Free settling / < /1 0:1 ÿ 0:3 g=l ws
1gD2s
s ÿ 18m
;
4
where m is the kinematic viscosity for clear water, s the speci®c gravity of suspended sediment and Ds the ¯oc diameter. (ii) Flocculation settling /1 < / < /2 0:3 ÿ 10 g=l ws k1 /4=3 ;
5
where k1 is an empirical coecient. (iii) Hindered settling / > /2 4:66
ws ws0 1 ÿ k2 / ÿ /2
;
6
where ws0 is the value of ws at the concentration level /2 , and k2 the inverse of /s ÿ /2 , where /s is the concentration at which ws 0. In estuarine and coastal waters, the horizontal length scale is generally much larger than the vertical scale, and hence an operator splitting algorithm may be used to separate the three-dimensional, advective±diusion equation into a horizontal, two-dimensional equation and a vertical, one-dimensional equation. Firstly, the two-dimensional advective±diusion equation was solved horizontally, and secondly, the one-dimensional vertical advective±diusion equation was solved vertically. To be consistent with the three-dimensional layerintegrated hydrodynamic model, the horizontal twodimensional equation was also integrated through the layers to give a layer-integrated two-dimensional equation. The two-dimensional QUICKEST scheme extended from LeonardÕs (1979) one-dimensional QUICKEST scheme was used to solve the horizontal two-dimensional layer-integrated equation. The modi®ed one-dimensional ULTIMATE algorithm was also used to prevent non-physical numerical oscillations, which frequently arise due to numerical dispersion when large concentration gradients exist (see Wu and Falconer, 184
1998). Another feature of the present model is that an additional source term associated with layer-integration of the free-surface ¯ow was introduced into the model, which was found to be of vital importance for mass conservation of pollutants (see Wu and Falconer, 1998). Boundary conditions The open boundary conditions were speci®ed by prescribing concentration data from ®eld measurements when in¯ow conditions prevailed, and with the concentration pro®les along the boundary being obtained by extrapolation using a ®rst-order upwind dierence scheme for out¯ow conditions. At bank (or closed) boundaries, normal derivatives of the concentration were set to zero. The one-dimensional, advective±diusion equation was solved using a non-uniform grid spacing in the vertical direction (see Fig. 1). Since the diusion term was very important for the vertical concentration distribution, and since some of the grid sizes were very small near the sea bed and water surface, the implicit power-law scheme of Patanker (1988) was adapted to avoid the use of very small time steps. For the sediment deposition and erosion rates, unlike two-dimensional depth-integrated models where these processes are modelled using sink and source terms, these ¯uxes are represented in the three-dimensional model using the following bed conditions: ÿ ws / ÿ ez
o/ qdep oz
when sb 6 sd
deposition; 7
ÿ ws / ÿ ez
o/ qero oz
when sb P se
erosion;
ÿ ws / ÿ ez
o/ 0 oz
when sd < sb < se
8
equilibrium; 9
where sb is the bed shear stress, sd the critical shear stress beyond which there is no further deposition, se the
Volume 37/Numbers 3±7/March±July 1998
critical shear stress for erosion, and qdep and qero represent the deposition and erosion rates, respectively, at the bed. In recent years, considerable eort has been made to study the mechanisms of cohesive sediment transport, and many experimental and ®eld studies have been carried out to investigate deposition and erosion rates of cohesive sediments (e.g. Krone, 1962; Mehta, 1973; Thorn and Parsons, 1977; Mehta, 1988; Parchure and Mehta, 1985; Thorn, 1981; Lick et al., 1995). The deposition rate proposed by Krone (1962) was used in this study: sb ; 10 qdep ÿws /b 1 ÿ sd where /b is the near-bed cohesive sediment concentration, with typical values of the critical shear stress sb for 2 deposition being 0.04±0.15 N=m (Krone, 1962; Mehta, 1973). Likewise, the erosion rate for soft natural mud can be represented by the following empirical expression (Raudkivi, 1990): qero 1=2 exp a sb ÿ se ; 11 qf where the empirical coecients are qf 0:0000042 and a 8:3 as suggested by Thorn and Parsons (1977). Typical values of the critical shear stress se for the ero2 sion of soft mud are 0.07±0.17 N=m . To complete the solution of the vertical, one-dimensional cohesive sediment transport equation, the following condition was applied at the free surface: w ÿ ws / ÿ ez
o/ 0 oz
12
which describes the condition that there can be no ¯ux across the free surface.
Humber Estuary Study General The Humber Estuary is one of the main estuaries in the UK, situated along the coast of northern England, and strategically important in terms of its ¯ow and solute input to the North Sea and its shipping links to mainland Europe. The main part of the estuary is about 62 km long, from Trent Falls to Spurn Head (see Fig. 2), and the tidal in¯uence extends for another 62 km up the River Ouse and about 72 km up the River Trent. As part of this study within the LOIS project, a comprehensive ®eld observation programme was undertaken and the model outlined herein was set up and applied to simulate cohesive sediment transport processes in the estuary from Trent Falls to just beyond Spurn Head. Field observations Observations of water elevations, currents, salinity, temperature and suspended particulate matter concentrations were made in the Humber mouth and plume
region of freshwater in¯uence during June 1995. Four tidal-cycle anchor stations were worked from R.V. Challenger, three of which were in the plume region: Stations SG13, SG24 and SG10. Station SG13 (53° 35.190 N, 0° 13.900 E) was occupied on 3 June between 05:17 and 18:45 h. The tidal range at Immingham was 4.9 m (mean tides). Station SG24 (53° 32.060 N, 0° 15.290 E) was worked on 4 June between 05:47 and 18:47 h. The tidal range was 4.5 m (large neaps). SG10 (53° 36.820 N, 0° 11.350 E) was worked on 7 June between 06:05 and 19:02 h. The tidal range was 3.8 m (neaps). The bed sediments at these stations were broken shells, coarse sand and larger material. The fourth Station SG23 (53° 35.710 N, 0° 2.170 E) was positioned within the estuary at Hawke Anchorage and was worked on 8 June between 07:21 and 20:53 h. The tidal range was 3.9 m (neaps). The bed sediments at this site were ®ne sand and smaller material. A tidal-cycle anchor station was worked from R.V. Sea Vigil, in the Trent Falls region of the upper Humber (53° 42.80 N, 0° 37.80 W). The station was occupied on 18 July between 09:04 and 21:54 h. The tidal range at Immingham was 5.3 m (small springs). The bed sediments were very ®ne sand and mud. The experimental procedure comprised water-column pro®ling from surface to bed at 1 or 2 m depth intervals, dependent upon the total water depth, and at approximately 0.5 h time intervals during each tidal cycle. The ship remained at anchor throughout the observations. A direct reading current meter was used to measure the current speed and direction. Direct reading temperature, conductivity, salinity and depth sensors were deployed, together with multi-range transmissometers for the determination of suspended particulate matter concentrations. Instruments were calibrated by comparison with data from laboratory analyses of water sampled during the measurements. Measurements in the plume area, just outside of the estuary and the model boundary, showed that suspended particulate matter, SPM, concentrations were small. At SG13, for example, the mean and standard deviation (SD) of the SPM for the whole tide was 13 4 mg lÿ1 . The mean and SD of the salinity was 33.8 0.7. There was a strong inverse relationship between SPM and salinity that indicated both that the dominant source of SPM was the estuary and that SPM transport was essentially conservative on the intra-tidal time-scale. Similar results were obtained at SG24 and SG10. The situation was dierent within the estuary at SG23. Here the SPM concentrations were signi®cantly greater and there was evidence of resuspension of sediments at the strongest current speeds. For the whole tide, SPM and salinity were 46 25 mg lÿ1 and 31.4 1.3, respectively. Measurements at Trent Falls showed that SPM concentrations were much greater in the upper reaches. For the whole tide, the means and SD of SPM and salinity were 2.4 2.6 g lÿ1 and 14.5 4.2, respectively. There was a signi®cant inverse correlation with salinity that 185
Marine Pollution Bulletin
Fig. 2 Map of the Humber Estuary and the locations of ®eld observation stations.
indicated that much greater SPM concentrations existed within the Ouse and Trent sub-estuaries. These turbidity maximum regions of extremely high SPM concentrations comprised suspended ®ne sediment (<63 lm), in which the majority of sediment existed as ¯ocs of clay and very ®ne silt particles. Model application and results The model area was represented horizontally using a mesh of 118 56 uniform grid squares, with a grid size of 500 m. Vertically, eight layers were used, with the thickness of the top layer being 4 m at mean water level, and with the other layers each being 3 m thick. The streamline boundary was located such that the boundary was as parallel as possible to the dividing streamline that separated the north±south alongshore current and the estuarine in¯ow±out¯ow. The water elevation recorded at Station SG13 was chosen as the seaward boundary condition to drive the tidal current. Due to the lack of volume ¯ux data at Trent Falls (Fig. 2), a water elevation boundary condition was prescribed at this landward boundary. A bed roughness height of 10 mm was assumed, as suggested by model calibration for a previous study (see Lin and Falconer, 1997). To achieve better predictions of the cohesive sediment distributions, accurate modelling of the hydrodynamic conditions in the estuary was extremely important. Since measurements taken at Trent Falls and Station SG13 were not simultaneous, i.e. the measurements at Station SG13 were taken for mean tide, whereas conditions at Trent Falls were taken during the spring tide, the water elevations at Trent Falls had to be adjusted to mean tide conditions. Calibration was carried out during the adjustment according to the water elevations and mean velocities measured at Hawke Anchorage. A comparison of the velocities for the dierent layers is given in Fig. 3. As can be seen from the results, the 186
general variation in time in the measured and predicted velocities at Hawke Anchorage shows a reasonable degree of agreement with regard to the depth-mean velocities. However, the variation over depth of the measured horizontal velocity components is generally much greater than the variation predicted from the numerical model. This under-prediction of the variation in the horizontal velocity ®eld is primarily thought to be due to the inclusion of a relatively simple turbulence model (i.e. the two-layer mixing length model) and the relative coarseness of the vertical grid, particularly near the bed boundary. Furthermore, the vertical advectionmomentum terms are represented using only a quasisecond order accurate scheme. These simpli®cations are currently being investigated further through another research programme. Suspended sediment concentrations are relatively low at Hawke Anchorage so that the free settling criterion is applicable [Eq. (4)]. Due to a lack of information on the primary particle sediment sizes at Hawke Anchorage, a sediment-speci®c gravity of 2.65 and a grain (micro¯oc) size of 20 lm was used in Eq. (4). Larger aggregates may be present, but these will have a lower speci®c gravity and settle relatively slowly. The critical shear stress for 2 deposition was set to: sd 0:07 N=m , and the critical 2 shear stress for erosion assumed to be: se 0:15 N=m . These values were in the ranges proposed respectively by Krone (1962) and Thorn and Parsons (1977), and were re®ned by trial and error. Fig. 4 shows a comparison of the ®eld-measured data and the model-predicted suspended cohesive sediment concentration variations with time. As the concentration variations with depth from the model predictions were small, only the depth mean concentration is shown in Fig. 4. The results generally show reasonable agreement between the model predictions and the ®eld-measured data in the regime of deposition. However, in the regime of erosion, the ®eld data showed a relatively large
Volume 37/Numbers 3±7/March±July 1998
Fig. 3 Comparison of predicted and ®eld measured velocities at Hawke anchorage, taken on 8 June 1995.
Fig. 4 Comparison of predicted and ®eld measured cohesive sediment concentrations at Hawke anchorage, taken on 8 June 1995.
variation with depth that has not been reproduced in the model predictions. This discrepancy is primarily thought to be due to the exclusion of non-cohesive sediment transport predictions in the model simulations reported herein. For the case of cohesive sediment transport predictions, once the erosive shear stress has been exceeded then the ®ne particulate matter will rapidly be distributed fairly uniformly through the water column in a strongly mixed estuary, because the settling velocity is very small. In contrast, for the case of non-cohesive sediment transport predictions, the settling velocity is much greater and the vertical variation much more pronounced. For example, previous simulations of the vertical variation of non-cohesive suspended sediment ¯uxes in the Humber Estuary showed very good agreement with ®eld data through the water column at Sunk Channel and Middle Shoal (Fig. 2; see Lin and Falconer, 1996). Hence, the disparity between the measured
and predicted sediment concentration distributions through the water column is thought to be due to the existence of a mixture of cohesive and non-cohesive sediment particles occurring in-situ during the erosion phase, whereas only cohesive sediment predictions have been included in the numerical model. This degree of variation between the predicted and measured suspended sediment concentration distributions is currently being investigated further in a similar model study for the Mersey Estuary, where the modelling of sediment mixtures is being considered.
Discussion and Conclusions A re®ned, three-dimensional, layer-integrated numerical model has been applied to simulate cohesive sediment transport processes in the Humber Estuary, UK. A combined layer-integrated and depth-integrated 187
Marine Pollution Bulletin
scheme was used to solve the hydrodynamic equations. The original three-dimensional advective±diusion equation was split into a two-dimensional horizontal and a one-dimensional vertical equation to implement the variation between the horizontal and vertical scales. A split algorithm was used to separate the three-dimensional transport equation into a horizontal, twodimensional equation and a vertical, one-dimensional equation due to the dierent length scales in the horizontal and vertical planes. To be consistent with the three-dimensional layer-integrated hydrodynamic model, the two-dimensional horizontal advective±diusion equation was integrated over the layers to give the layerintegrated advective±diusion equation. A conservative form of the horizontal, two-dimensional, layer-integrated advective±diusion equation was used to ensure mass conservation of the cohesive sediment ¯uxes. A two-dimensional modi®ed QUICKEST scheme was developed and used in this study to achieve high numerical accuracy. To avoid numerical oscillations, a modi®ed one-dimensional universal limiter ULTIMATE scheme was also used. Comparisons between the measured and predicted velocities at Hawke Anchorage generally showed good agreement between the depth-mean velocities. However, the measured velocities showed a greater variation about the mean value through the water column. This was thought to be primarily due to the relatively simple turbulence model and the coarse grid resolution near the bed. The model is currently being extended to include a version of the two-equation k±e turbulence model in the vertical. Comparisons of the suspended cohesive sediment concentrations were also undertaken at Hawke Anchorage and showed good agreement for the deposition phase. However, during the erosion phase the depthmean concentrations were generally in good agreement, but again the variation over the water column was much greater for the measured data than for the predicted results. This disparity in the vertical concentration distribution during the erosion phase was primarily thought to be due to mixtures of cohesive and non-cohesive sediment ¯uxes occurring in-situ, whereas only cohesive sediment transport predictions were included in the model. These results have indicated that further research needs to be undertaken on predicting ¯uxes of sediment mixtures. This is LOIS publication number 551 of the LOIS Thematic Research Programme, carried out under a Special Topic Award from the Natural Environment Research Council. The authors are grateful to NERC for supporting this study. Ariathurai, R. and Krone, R. B. (1976) Finite element model for cohesive sediment transport. Journal of Hydraulic Engineering Division, ASCE 102, 323±338. Cole, P. and Miles, G. V. (1983) Two-dimensional model of mud transport. Journal of Hydraulic Engineering Division, ASCE 109, 1± 12.
188
Davies, A. M. (1983) Formulation of a linear three-dimensional hydrodynamic sea model using a Galerkin-eigenfunction method. International Journal for Numerical Methods in Fluids 3, 33±60. Falconer, R. A. and Chen, Y. (1991) An improved representation of ¯ooding and drying and wind stress eects in a 2-d numerical model. Proceedings of the Institution of Civil Engineers, Part 2, Research and Theory 91, 659±687. Falconer, R. A. (1993) An introduction to nearly horizontal ¯ows. In Coastal, Estuarial and Harbour Engineers' Reference Book, ed. M. B. Abbott and W. A. Price, pp. 27±36. E & FN Spon Ltd., London. Falconer, R. A. and Chen, Y. (1996) Modelling sediment transport and water quality processes on tidal ¯oodplains. In Floodplain Processes, eds. M. G. Anderson, D. E. Walling and P. D. Bates, pp. 361±398. Wiley, New York. Hayter, E. J. and Mehta, A. J. (1986) Modelling cohesive sediment transport in estuarine waters. Applied Mathematical Modelling 10, 294±303. Krone, R. B. (1962) Flume studies of the transport of sediment in estuarial processes. Final Report, Hydraulic Engineering Laboratory and Sanitary Engineering Research Laboratory. University of California, Berkeley. Leendertse, J. J., Alexander, R. C. and Liu, S. K. (1973) A threedimensional model for estuaries and coastal seas: Volume I. Principles of Computation. Report No. R-1417-OWRR, Rand Corporation, Santa Monica. Lick, W., Xu, Y. J. and McNeill, J. (1995) Resuspension properties of sediments from the Fox, Sagina and Bualo rivers. Journal of Great Lakes Research 21, 257±274. Lin, B. and Falconer, R. A. (1996) Numerical modelling of threedimensional suspended sediment for estuarine and coastal waters. Journal of Hydraulic Research, IAHR 34, 435±456. Lin, B. and Falconer, R. A. (1997) Three-dimensional layer-integrated modelling of estuarine ¯ows with ¯ooding and drying. Estuarine, Coastal and Shelf Science 44, 737±751. Lin, P., Huan, J. and Li, X. (1983) Unsteady transport of Suspended load at small Concentrations. Journal of Hydraulic Engineering, ASCE 109, 86±98. Markofsky, M., Lang, G. and Schubert, R. (1986) Suspended sediment transport in rivers and estuaries. In Lecture Notes on Coastal Engineering Studies, Physics of Shallow Estuaries and Bays, ed. J. Van de Kreeke, pp. 210±227. Mehta, A. J. (1973) Deposition behaviour of cohesive sediment, Ph.D. Thesis, University of Florida, Gainesville. Mehta, A. J. (1988) Laboratory studies on cohesive sediment deposition and erosion. In Physical Processes in Estuaries, ed. J. Dronkers and W. van Leussen, pp. 427±445. Springer, New York. Mehta, A. J. (1993) Hydraulic behaviour of ®ne sediment. In Coastal, Estuarial and Harbour Engineers' Reference Book, ed. M. B. Abbott and W. A. Price, pp. 577±584. E & FN Spon Ltd., London. Nicholson, J. and OÕConnor, B. A. (1986) Cohesive sediment transport model. Journal of Hydraulic Engineering, ASCE 112, 621±640. Odd, N. V. M. and Owen, M. W. (1972) A two-layer model for mud transport in the Thames Estuary. Proceedings of the Institution of Civil Engineers, London, Supplement (ix), pp. 175±205. Parchure, T. M. and Mehta, A. J. (1985) Erosion of soft cohesive sediment deposits. Journal of Hydraulic Engineering, ASCE 111, 1308±1326. Raudkivi, A. J. (1990) Loose Boundary Hydraulics. 3rd edn. Pergamon Press, Oxford. Rodi, W. (1984) Turbulence Models and their Application in Hydraulics. 2nd edn. IAHR Publication, Delft. Scarlatos, P. D. (1981) On the numerical modelling of cohesive sediment transport. Journal of Hydraulic Research, IAHR 19, 61± 68. Shen, Y. P. (1987) On modelling three-dimensional estuarine and marine hydrodynamics. In Three-Dimensional Models of Marine and Estuarine Dynamics, ed. J. C. J. Nihoul and B. M. Jamart, pp. 35± 44. Elsevier Oceanography Series, Elsevier, Amsterdam. Thorn, M. F. C. and Parsons, J. G. (1977) Properties of Grangemouth Mud. Hydraulic Research Station, Rep. No. EX 781, Wallingford, UK. Thorn, M. F. C. (1981) Physical processes of siltation in channels. Proceedings of the Conference on Hydraulic Modelling Applied to Maritime Engineering Problems, pp. 47±55. Institution of Civil Engineers, London. Wu, Y. and Falconer, R. A. (1998) Two-dimensional ULTIMATE QUICKEST scheme for pollutant transport in free-surface ¯ows.
Volume 37/Numbers 3±7/March±July 1998 Third International Conference on Hydroscience and Engineering. Cottbus/Berlin, Germany.
Yost, S. A. and Katopodes, N. D. (1996) 3D numerical modeling of cohesive sediments. Proceedings of the International Conference on Estuarine and Coastal Modeling, pp. 478±489.
189