Modelling of hydrodynamics and cohesive sediment transport in
Tanshui River estuarine system, Taiwan
Wen-Cheng Liu
a,*, Ming-Hsi Hsu
b, Albert Y. Kuo
baHydrotech Research Institute, National Taiwan University, No. 158, Chou-Shan Road, Taipei 10617, Taiwan bDepartment of Bioenvironmental Systems Engineering, National Taiwan University, Taipei 10617, Taiwan
Abstract
Alaterally averaged two-dimensional numerical model is used to simulate hydrodynamics and cohesive sediment transport in the Tanshui River estuarine system. The model handles tributaries as well as the main stem of the estuarine system. Observed time series of salinity data and tidally averaged salinity distributions have been compared with model results to calibrate the turbulent diffusion coefficients. The overall model verification is achieved with comparisons of residual currents and salinity distribution. The model reproduces the prototype water surface elevation, currents and salinity distributions. Comparisons of the suspended cohesive sediment concentrations calculated by the numerical model and the field data at various stations show good agreement.
The validated model is applied to investigate the tidally averaged salinity distributions, residual circulation and suspended sediment concentration under low flow conditions in the Tanshui River estuarine system. The model results show that the limit of salt intrusion in the mainstem estuary is located at Hsin-Hai bridge in Tahan Stream, 26 km from the River mouth under Q75flow.
The null point is located at the head of salt intrusion, using 1 ppt isohaline as an indicator. The tidally averaged sediment con-centration distribution exhibits a local maximum around the null point.
Ó 2002 Elsevier Science Ltd. All rights reserved.
Keywords: Tanshui River estuary; Estuarine hydrodynamic; Cohesive sediment transport; Salinity distribution; Residual circulation; Calibration and verification
1. Introduction
Estuaries, the transition zones between the rivers and the sea, show a complex dynamic behaviour. Salt water intrusion and mixing of salt and fresh water can led to a complex stratification pattern with associated density currents. Cohesive sediment transport in estuaries can led to the occurrence of a turbidity maximum, which induces problems, both ecological and economical. The movement of cohesive sediments can cause siltation to occur in waterways and harbor docks, or erosion of estuarine banks. To deal with these problems, expensive dredging operations or bank protections are generally needed. The contaminated materials from industrial and municipal effluents, or accidental oil spills, can also re-lease heavy metals, mineral oils and other toxic con-taminants, which may then be absorbed onto fine
sediments and be available for resuspension by strong tidal currents, short wave action and dredging opera-tions. This desorption of contaminants from their par-ticulate phase can impact significantly on the ecological balance of shallow, near-shore areas. Therefore, accu-rate predictions of cohesive sediment transport pro-cesses in estuarine and coastal waters are of vital importance to environmental management.
Hydrodynamic processes are the most important mechanisms involved in sediment transport. They ad-vect the suspended sediments, provide the force need to erode the bed and, through turbulence, play a major role in the flocculation of cohesive sediments. Relatively large velocities generally occur in tidal estuaries. The sediments are very mobile in estuaries. They are eroded and transported upriver during flood, deposited during slack water, eroded again and transported downriver during ebb and redeposited during next slack water, to restart their movement in the forthcoming tidal cycle.
In recent years, the advent of computer has opened up the possibility of modeling cohesive sediment transport
*Corresponding author.
E-mail address:wcliu@hy.ntu.edu.tw(W.-C. Liu).
0025-326X/02/$ - see front matterÓ 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 5 - 3 2 6 X ( 0 2 ) 0 0 1 6 0 - 1
numerically. Odd and Owen (1972) developed a one-dimensional two-layer model considering a deeper thin lower layer with constant depth and an upper layer of varying thickness. Erosion and deposition rates were based on the formulations proposed by Krone (1962) and Partheniades (1965). One-dimensional mod-els have been frequently used to simulate sediment transport and large-scale morphological changes in river (De Vries et al., 1989), in tidal channels (Dry and Evans, 1989) and for the simulation of lutocline formation in estuaries (Ross and Mehta, 1989; Smith and Kirby, 1989).
An early two-dimensional vertically integrated model was present by OÕConnor (1971). Ariathurai and Krone (1976) presented a finite element model adopting trian-gular elements with a quadratic approximation for the concentration and a Galerkian weighted residual method. Mulder and Udink (1991) presented a two-di-mensional model for the Western Scheldt estuary (The Netherlands) in which the tide and wind wave models were combined to produce stationary wave fields. Barros (1996) presented a two-dimensional depth-averaged formulation of the transport equation. The model is composed of two separate, though coupled, modules: one to describe water column dynamics (advection– dispersion, erosion, deposition, flocculation), and an-other to describe bottom consolidation mechanisms. The numerical implementation of the model consists of a finite-element solution of the transport equation on an Eulerian–Lagrangian framework.
Kuo et al. (1978) developed a laterally averaged two-dimensional model to investigate the observed turbidity maximum in the Chesapeake Bay tributary estuaries. They demonstrated that the convergent bottom residual flow near the head of salt intrusion may create the ob-served maximum there. Li et al. (1994) developed a coupled width-integrated hydrodynamic and sediment transport model for the Gironde estuary (France). A turbulent-closure model is described to compute the turbulent viscosity and diffusion coefficients. Bottom exchange was calculated based on the classical formu-lations.
Although a fully three-dimensional numerical model is generally desirable from a physics point of view, the verification of sediment transport models in estuaries and adjacent shallow embayments has been often found to be more reliable for the case of two-dimensional models as opposed to their three-dimensional counter-parts (Van Rijn et al., 1990). This is due to high the cost and the technological limitations inherent to the acqui-sition of data with the quality and quantity needed to validate a three-dimensional model. Field measurements of sediment transport, erosion and deposition in a three-dimensional framework are rare in the literature, and when available they have been carried out for short periods and within limited areas of an estuary.
In this paper a two-dimensional, laterally averaged finite difference model of hydrodynamics and cohesive sediment is used to determine the respective role of several processes on the sediment transport character-istics in the Tanshui River estuarine system. Aconser-vative form of the two-dimensional advective-diffusion equation has been used to ensure conservation of the cohesive sediment mass. The model has been set up and applied to predict cohesive sediment transport in the Tanshui River estuarine system. The model predictions and field measured data were compared to evaluate the model reliability.
2. The Tanshui River estuarine system
The Tanshui River is formed by the confluence of the Tahan Stream, Hsintien Stream and Keelung River (Fig. 1). The downstream portions of all three tributaries are influenced by tide, and subjected to sea water intrusion. Together, they form the largest estuarine system in Taiwan, with its drainage basin including the capital city of Taipei. Its drainage area encompasses 2728 km2 with a total channel length of 327.6 km. The major portion of the estuarine system, upriver of Kuan-Du (Fig. 1), lies within the Taipei basin, while the stretch downriver is confined by high mountains on both sides. Because of the presence of mountains and narrowness of the river, there is no significant wind-induced current. Except the occasional storm surges induced by hurricanes, the major forcing mechanisms of the barotropic flows are astronomical tide at the river mouth and river dis-charges at upriver ends. Semi-diurnal tides are the principal tidal constituents, with a mean tidal range of 2.22 m and a spring tidal range of 3.1 m. The average river discharges are 62.1, 72.7, and 26.1 m3/s, respec-tively, in the Tahan Stream, Hsintien Stream and Kee-lung River. In addition to the barotropic flows forced by tide and river discharges, the baroclinic flow forced by sea water intrusion is another important transport mechanism in the Tanshui River estuarine system. The mean annual river flow and the total solid matters dis-charged from the upland basin have been estimated to be 7044 106 m3/y and 11:45 106 t/y, respectively, over the last 40 years (Hydrological Year Book of Tai-wan, 1996). Most of the solid matters are non-cohesive sediments and they reach the tidal reach only during the times of flood events.
3. The mathematical model 3.1. Approach
The tidal-time model consists of a systematic se-quence of numerical procedures designed to describe
and to simulate the main movement of water and sedi-ment in an estuarine system. The chief variables of nu-merical computation are tidal height, current, salinity, and suspended sediment concentration in response to the basic external factors: tidal wave propagation, river flow, and river-borne sediment fluxes. Inasmuch as the forces and variables are continuously interacting and changing with time, the model is time dependent.
Since the Tanshui River estuarine system is long and narrow and suspended sediment concentrations change with distance seaward and with depth, the model is two-dimensional. It is assumed that all properties are uni-formally distributed laterally across the estuary at their respective average values.
To calculate the movement of water flow and suspended sediment distribution, the time dependent, Fig. 1. (a) The map of the Tanshui River estuary, and (b) the model segments. (Lines cutting across the river are model transects at 1.0 km interval; numbers are model segment numbers.)
three-dimensional equations of motion and the equation of mass balance are first integrated laterally to obtain two-dimensional equations. Next, the two-dimensional equations are integrated with depth over the height of successive layers, and the resulting equations for each layer are written in finite difference forms. The finite difference equations are then solved numerically with prescribed initial and boundary conditions.
The numerical computation yields the time-varying tidal height, current, salinity, and suspended sediment concentration distributions throughout the estuary. Thus, these time variations display the response of sa-linity and sediment distributions to tidal current fluc-tuations such as the resuspension of sediment particles during increasing flood or ebb current strength as well as deposition during slack current periods. Sediment flux at any transect may be calculated by multiplying the hor-izontal velocity with sediment concentration and inte-grating over the cross-section. Residual values can be obtained for each variable by averaging the respective variable over a complete tidal cycle. The average sedi-ment concentration field defines the turbidity maximum if the hydrographic conditions lead to the development of a maximum.
3.2. Hydrodynamic model
The hydrodynamic model used to predict the water surface elevations, longitudinal and vertical velocities, and salinity involves solving the governing equations of fluid flow. Generally, the flow was assumed to be iso-thermal and the vertical acceleration was assumed to be small compared to the gravitational acceleration, yield-ing a hydrostatic pressure distribution. Tidal forcyield-ing is prescribed at the downstream boundary, whereas the flow rates are imposed at the tidal limits at upstream boundaries. In addition to the momentum and conti-nuity equations, the model also solves transport equa-tion for salt and an equaequa-tion of state to include the baroclinic effects. The turbulent shear stresses are 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, respectively. Full details of the governing equations used in the model are given in Hsu et al. (1999).
The model solves the continuity, momentum and mass-balance equations by a two-time level, finite dif-ference scheme with spatially staggered grid. The im-plicit treatment of the vertical diffusion terms results in a tri-diagonal matrix in the vertical direction. The QUICKEST scheme is used to express the finite-differ-ence form of the advective term. It is based on a con-servative control volume formulation and estimates the cell wall concentration with a quadratic interpolation using concentrations in two adjacent cells and that at the
next upstream cell. The method of the solution is de-tailed in Park and Kuo (1993).
3.3. Sediment transport model
The following advection–diffusion, partial differential equation is used to describe the suspended sediment transport processes: oðBCÞ ot þ oðuBCÞ ox þ o½ðw WsÞBC oz ¼ o ox KxB oC ox þo oz KzB oC ox Bd þ Br þ qCt ð1Þ where t is time, x and z are longitudinal and vertical coordinates, respectively. C is the cohesive sediment concentration, B is the river width as a function of x and z, u and w are the velocity components in the x and z directions, respectively, Ws is the apparent sediment
settling velocity, d is sediment deposition rate, r is sed-iment resuspension rate, q is tributary inflow through unit area of x–z plane, Ct is the sediment concentration
of tributary flow, Kx and Kz are turbulent diffusion
co-efficients in the x and z directions, respectively. The deposition rate, d, and the resuspension rate, r, are zero except at the bottom interface.
The Munk–Anderson-type formulations is used to calculate the vertical turbulent diffusion coefficient, Kz. The model used constant values for Kxand it is adjusted, within the range of 104 to 106 cm2/s, through model calibration (Hsu et al., 1999).
3.4. The settling velocity
The vertical transport is due to vertical advection, particle settling and turbulent diffusion. The hydrody-namic model computes the vertical velocity and the turbulent diffusivity. The settling velocity depends on gravitational force, and on the vertical shear due to settling movement. The gravitational force depends on the density of each individual particle (terrigenous or biological) forming the flocs, and on the floc porosity occupied by water.
As a result of floc aggregation due to inter-particle collisions and the surface electro-chemical forces, cohe-sive sediments settle by flocs rather than by individual particle. The frictional forces depend on the form of the floc and on the Reynolds number of the surrounding flow during settling. For very small particles in laminar flow the ratio between the gravitational and the fric-tional forces is proporfric-tional to the square of the particle diameter of the floc. So, the settling velocity is expected to increase with floc size.
Sverdrup et al. (1942) indicated that StokesÕ law can be considered valid for spheres of diameters smaller than
62.5 lm in quiescent waters. As pointed out by Sverd-rup, turbulent diffusion in the water column can cause an increase or a decrease in the downward transport due to settling. Since StokesÕ law was derived without con-sideration of turbulence, Ws should be altered to reflect estuarine turbulent conditions. Unfortunately, experi-ments in the Tanshui estuary involving small size par-ticles (<1 mm) have not been conducted to discover the relationships between turbulence and settling velocities predicted by StokesÕ law. Thus, the only theoretical basis on which settling velocity for spheres of diameters smaller than 62.6 lm may be calculated is StokesÕ law: Ws¼ 1 18m qs q q gD 2 ð2Þ
where qsis the density of sediment particle, q is density
of water, m is the kinematic viscosity of water, D is StokesÕ diameters of particles.
StokesÕ law shows that the settling velocity is not constant for all sediment particles. Assuming qs is con-stant for particles of all sizes and neglecting the density variation of estuarine water, the settling velocity may be considered to be proportional to the square of particle diameter, i.e. WsðDÞ ¼ KD2 ð3Þ where K¼ 1 18m qs q q g
To evaluate the vertical sediment fluxes, terms con-taining WsC, the size dependence of settling velocity has to be taken into account. Let cðDÞ be the concentration of the sediment particles with a StokesÕ diameter in the range D and Dþ dD, then
WsC¼ Z 1 0 WsðDÞcðDÞ dD ¼ K Z 1 0 D2cðDÞ dD ð4Þ Statistics give the variance of diameters of sediment particles as h2¼ R1 0 ðD DmÞ 2 cðDÞ dD C ð5Þ
where the total concentration C¼
Z 1 0
cðDÞ dD and the mean diameter Dm¼ 1 C Z 1 0 DcðDÞ dD then Eq. (5) gives Z 1 0 D2cðDÞ dD ¼ Cðh2þ D2 mÞ ð6Þ and thus WsC¼ KðD2mþ h 2ÞC ð7Þ
which is used in the mass-balance equation for the set-tling velocity. The variation of K with time is small, and therefore, K is assumed to be constant. The same as-sumption is made for h2 and D2
m. Sampling along the
extent of the Tanshui River estuary on February 2 and 23, 1999 results in an average diameter of 16.5 lm and a mean variance of 78.8 lm2(Figs. 2 and 3). The average
settling velocity of 0.026 cm/s is adopted in the model simulation.
3.5. The erosion model
Erodibility of a cohesive bed is driven by shear, but also depends on bottom cohesive nature, which in turn depends, in a poorly understood way, on clay mineral-ogy and on the geochemistry and miscrobiological processes occurring in the bottom. Some authors argued that it should also depend on the salinity (Hayter and Mehta, 1986). However, no dependency laws have yet been advanced.
Fig. 2. Mean particle size of suspended sediment in the lower estuary (each data point represents mean of one sample).
Fig. 3. Variance of particle size of suspended sediment in the lower estuary (each data point represents variance of one sample).
Again, a useful correlation must depend only on the variables calculated by the model and on the parameters which may be determined by experiments or by model calibration. The erosion algorithm used in this work is based on the classical approach of Partheniades (1965). Erosion occurs when the bed shear stress exceeds the threshold of erosion. The flux of eroded matter is given by r¼ M s sE 1 for s > sE ð8aÞ r¼ 0 for s 6 sE ð8bÞ
where s is the bed shear stress, sE is the critical shear
stress for erosion and M is the erosion constant. The parameter M depends on the physico-chemical characteristics of bottom sediment. In the model simu-lation of Tanshui River estuary, a value of 3 106 kg m2s1 is used. As a general rule, bottom sediments are mixtures of cohesive and non-cohesive sediments; this parameter must also account for that, and a gra-dient must be expected in the estuary.
Critical shear stress for erosion is a function of the degree of compaction of bottom sediments measured by the dry density of the bottom sediments: ratio between the mass of sediment (after extraction of the interstitial water at 105 °C) and its initial volume. Delo (1988) proposed the formulation which is used in this study sE ¼ A1ðqdÞ
E1 ð9Þ
where qdis the dry density of bed sediments, and A1and
E1 are coefficients depending on mud type.
Using the bed sediments collected from the Tanshui estuary, experiments were conducted to understand the mechanism of sediment erosion. Based on the definition proposed by McNeil et al. (1996), the critical conditions of erosion were investigated quantitatively. For each bed sediment sample, the critical flow velocity, critical shear stress, and the dry density were measured. Figs. 4 and 5
present the experimental results with samples collected near Tanshui River mouth and Zhu-Wuei station. The equations of sE¼ 0:83q1:04d and sE¼ 0:98q1:14d were de-rived by liner regression of data from samples near river mouth and Zhu-Wuei station, respectively. The results show that the dry density representing the condition of consolidation is highly correlated to the critical velocity as well as to the critical shear stress, and an exponential function may be established for model applications. From the experiment, we found that dry density of un-processed samples range from 0.15 to 0.2 g/cm3 in the Tanshui estuary. The corresponding sE of 0.1 N/m2 is adopted in the model simulation.
3.6. The deposition model
The deposition algorithm, like the erosion algorithm, is based on the assumption that deposition and ero-sion never occur simultaneously. An algorithm was first proposed by Krone (1962) and later modified by Odd and Owen (1972). The algorithm is based on the as-sumption that a particle reaching the bottom has a probability of remaining there that varies between 0 and 1 as the bottom shear stress varies between its upper limit for deposition and zero, respectively. The deposi-tion rate is calculated as the product of the settling flux and the probability of a particle to remain on the bed d ¼ CBWs 1 s sD for s < sD ð10aÞ d ¼ 0 for s P sD ð10bÞ
where sDis the critical shear stress for deposition and CB
is the sediment concentration at the bottom layer. The critical shear stress for deposition, sD, depends mainly on the size of the floc. Bigger flocs have higher probability of remaining on the bed than smaller flocs. Nevertheless, previous work suggested that a constant Fig. 4. Relationship between critical velocity for erosion and sediment
dry density.
Fig. 5. Relationship between critical shear stress for erosion and sed-iment dry density.
value is a reasonable approximation. The value of 0.05 N/m2 is used in the model simulation.
4. Model calibration and verification 4.1. Hydrodynamics
The hydrodynamic model is supplied with data de-scribing the geometry of the Tanshui River system. The geometry in the vertical two-dimensional model is rep-resented by width at each layer at the center of each grid cell. Afield survey in 1994 and 1995 by Taiwan Water Conservancy Agency collected the cross-sectional pro-files at about 0.5 km interval along the tidal portion of the river. These profiles were used to schematize the river system. The three branches are divided into 33, 14, and 37 segments, respectively, with a uniform segment length of 1.0 km. The vertical layer thickness is 1 m for all layers, except the surface layer which is variable, with 2 m at mean sea level. The maximum number of layers is 10 at the deepest section of the river. To avoid upstream boundary conditions being influenced by tide, the computational domain of a model should be extended to or beyond the tidal limits in the tributaries as well as the mainstem. For a tidal river with a well-defined fall line, the tidal limit is located at fixed point. It is the ideal location as upstream boundary of model domain. The boundary condition there may be specified with fresh-water discharge as mass input without contribution to momentum. For a tidal river with bottom elevation rising gradually above mean sea level, the tidal limit moves up or down the river in response to a decrease or increase of freshwater discharge. The model limit should be at a location where the bottom elevation is higher than high tide level. The freshwater discharge contrib-utes both mass and momentum at this type of upstream boundary condition. The Tanshui River system is of the second type of boundary condition. The hydrodynamic model was also expanded to include the capability to simulate tributaries as well as the mainstem of estuaries (Hsu et al., 1996, 1997, 1999). In applying to the Tan-shui River system, the model treats the TanTan-shui River and Tahan Stream as the mainstem and Hsintien Stream and Keelung River as tributaries. Details of hydrody-namic model calibration and verification have been presented in Hsu et al. (1996, 1997, 1999). In the fol-lowing, a brief description is presented.
ManningÕs friction coefficient and the coefficients for turbulent diffusion terms are important calibration pa-rameters affecting the calculation of surface elevation, current velocity and salinity distribution. The prelimi-nary calibration of ManningÕs friction coefficient used a single constituent tide, M2, to reproduce the longitudinal distribution of mean tidal range. The results show that the average absolute values of differences (mean
abso-lute difference) between computed and measured tidal ranges is 3.4 cm, while the root-mean-square difference is 3.8 cm, both of which are less than 2% of mean tidal range. The fine-tune calibration compares the along-river variations of amplitude and phase of the individual tidal constituent, using a nine-constituent tide to force the system at the river mouth. The results also compare times of high tide and low tide between computation and observation. The differences in the model results and observations are less than 10 min. The calibrated model has the ManningÕs friction coefficient ranging from 0.032 to 0.026 in the Tanshui River–Tahan Stream, 0.015 for the Hsintien Stream, and from 0.023 to 0.016 for the Keelung River.
The Taiwan Water Conservancy Agency conducted three intensive field measurements on April 12, June 24, 1994 and April 14, 1995. Half-hourly current speed and hourly temperature and salinity were measured by personel on boats for a period of 13 daylight hours in each case. Velocity data were taken with handheld cur-rent meters that measured curcur-rent magnitude but not the direction. The data were recorded by hand, and ‘‘ebb or flood’’ direction was noted with visual observation. Except for the upriver transects close to tidal limits, the measurements were made at several depths in the ver-tical direction and several points in the transverse direction. To compare with the results of the vertical two-dimensional model, the data in the transverse di-rection were averaged to obtain the laterally averaged values. Fig. 6 presents the comparison of the time series data of the longitudinal velocity on April 14, 1995. It verifies that the flow is properly simulated by the model. The mean absolute difference and the root-mean-square difference between computed and measured data range from 8.0% to 35.0% of peak velocities at the corre-sponding stations.
The calibration of mixing processes and verification of flow field was conducted through model simulation of prototype conditions during the period March 15–Sep-tember 30, 1994. The model was driven by three time-varying boundary conditions: measured daily freshwater inflows through the upstream boundaries, hourly tidal elevation and linearly interpolated salinity from mea-sured data at the river mouth. The mixing processes are modeled with turbulent diffusion terms. The mixing length concept is used to calculate eddy viscosity and diffusion coefficient in the vertical direction. The for-mulations similar to those proposed by Pritchard (1960) are used in the model. The tidally averaged salinity and time series data measured on April 12 and June 24, 1994 were used to calibrate the constants in the formulation of vertical diffusion coefficient. The results show that the mean absolute difference and the root-mean-square difference between the computed and hourly measured data range from 0.13 to 1.76 ppt and from 0.16 to 2.45 ppt, respectively, at various stations. The computed
water surface elevations and current velocities were compared with time series data at various stations to verify the barotropic flow component. The measured and computed longitudinal velocity at selected stations were averaged over tidal cycle. The vertical profiles of measured and computed residual velocity were com-pared and, also verified with analytical models of Han-sen and Rattray (1965). The results have been reported by Hsu et al. (1999).
The modelÕs ability to predict mass transport was verified with a simulation of salinity distribution from March 15 to April 30, 1995. The model predictions agree well with the field measurements on April 14, 1995. Fig. 7 shows samples of the time series comparison at se-lected stations. The results show that the mean absolute difference and the root-mean-square difference between the computed and hourly measured data range from 0.12 to 1.44 ppt and from 0.13 to 1.62 ppt, respectively, at various stations. The calibrated and verified values of the longitudinal dispersion coefficients are Kx¼ 28 105 cm2/s for the Tanshui River–Tahan Stream, Kx¼ 3:5 105cm2/s for the Hsintien Stream, and Kx¼ 5 104 cm2/s Keelung River.
4.2. Cohesive sediment transport
Measurements of the suspended sediment concen-trations in the Tanshui River system during 1994 and
1995 were used to calibrate and verify the cohesive sediment transport model. The data were provided by Taiwan EPA, which include the measured suspended sediment concentrations at the upstream and down-stream boundaries. They were interpreted in time and used for model boundary conditions. The model was run for two one-year simulations. The daily averages of the model results were compared with the values from field observations at various stations. Figs. 8 and 9 present the comparisons for calibration and verification runs of 1994 and 1995, respectively. Since there is no record of the depth of field sampling, we can only judge the model performance by the fact that the predicted surface and bottom sediment concentrations generally bracket the field data. It is apparent that the field ob-servations missed most of the peak concentrations pre-dicted by the model. These prepre-dicted peak sediment concentrations are the results of storm induced high freshwater flow coming from upriver, while almost all of the data were collected under fair weather condition. Nevertheless, the results show reasonable agreement between the model predictions and field observations.
Figs. 8 and 9 provide a qualitative comparison of model predictions and field observations. This tradi-tional assessment of model adequancy, the perceived agreement between predictions and observations, de-pends upon the viewpoint and experience of assessors. In order to render the evaluation of the model less Fig. 6. Model verification: comparison of velocity time series at (a) Kuan-Du bridge, (b) Taipei bridge, (c) Hsinhai bridge, and (d) Pa-Ling bridge.
Fig. 7. Model verification: comparison of salinity time series at (a) Kuan-Du bridge, (b) Taipei bridge, (c) Hsinhai bridge, and (d) Pa-Ling bridge.
Fig. 8. Model calibration: comparison of predicted and field measured 1994 cohesive sediment concentrations at (a) Kuan-Du bridge, (b) Chrong-Yang bridge, (c) Pa-Ling bridge, and (d) Hwa-Jiang bridge.
subjective, quantitative assessment of model accuracy are desirable. No single measure or set of measures is universally accepted for this purpose. The mean abso-lute difference and the root-mean-square difference be-tween computed and measured data are often used to evaluate the model performance. They are presented in Table 1 to show the quantitative assessments of model calibration and verification. Scatterplots of point-by-point comparison of predictions and observations at various stations are also presented in Fig. 10. Asolid, diagonal line indicates the perfect one-to-one corre-spondence.
5. Model application
The calibrated and verified model was used to per-form model simulations in the Tanshui River estuarine
system to investigate the cohesive sediment transport. The values of all coefficients of the numerical model have been determined by the calibration and verification processes; no further adjustment to the coefficients was made. The time series data of surface elevation at the Tanshui River mouth, collected by the Taiwan Water Conservancy Agency, were examined and harmonic analyses were performed. Nine principal constituents were identified, they are M2 (12.42 h), S2(12 h), N2(12.9 h), K1 (23.93 hr), Sa(8765.32 h), O1 (25.82 h), K2 (11.97 h), P1(24.07 h), and M4(6.21 h). Amplitudes and phases of the constituents were specified to generate a time se-ries of surface elevation as the downstream boundary condition for a one-year (705 tidal cycles) model simu-lation run. The Q75discharges (the flow which is equaled or exceeded 75% of time) at the upstream boundaries of the three major tributaries were used. They are 8.15, 20.2, and 3.61 m3/s for the Tahan Stream, Hsintien
Fig. 9. Model verification: comparison of predicted and field measured 1995 cohesive sediment concentrations at (a) Kuan-Du bridge, (b) Chrong-Yang bridge, (c) Pa-Ling bridge, and (d) Hwa-Jiang bridge.
Table 1
Quantitative comparison of observed and modelled concentration of suspended sediments for calibration and verification of cohesive sediment model Stations Number of observations Average observed
concentrations (mg/l)
RMS differences (mg/l) Mean absolute differences (mg/l)
Kuan-Du bridge 46 39.8 5.3 4.4
Chrong-Yang bridge 46 28.4 4.3 3.6
Pa-Ling bridge 46 42.1 8.0 5.1
Stream, and Keelung River, respectively. Salinity at the Tanshui River mouth was set at 32 ppt under Q75 flow condition. In the sediment transport model, the char-acteristic concentration of 20 mg/l was used for the boundary condition at the river mouth. The upriver boundaries are located at the limits of the tide in the three tributaries, and a constant concentration of 40 mg/l were imposed. The model was run for one-year duration (705 tidal cycles).
Fig. 11 presents the average salinity distributions over two spring-neap cycles (i.e., 58 tidal cycles) in the Tan-shui River-Tahan Stream, Keelung River, and Hsintien Stream. The limits of salt intrusion are located at Hsin-Hai bridge in Tahan Stream, about 26 km from the Tanshui River mouth, 12 km from its mouth in the Keelung River and 4 km from its mouth in the Hsintien River. The extensive intrusion of saline water imposes a significant baroclinic forcing and induces a strong re-sidual circulatory system in the estuary. Fig. 12 shows the tidally averaged longitudinal velocity distributions (averaged over two spring-neap cycles) in the Tanshui River–Tahan Stream, Keelung River, and Hsintien Stream. The null points, where the upriver and down-river bottom flows converge, are evident at the heads of salt intrusion, i.e. 1 ppt isohaline, as an indicator. The downriver net velocity in the upper layer increases sea-ward despite the enlargement of river cross-section in
Fig. 10. Scatterplots for model calibration and verification. (The model predicted concentrations are average values over depth and tidal cycle.)
Fig. 11. Calculated salinity distributions averaged over 58 tidal cycles under Q75 flow condition in (a) Tanshui River–Tahan Stream, (b) Keelung River, and (c) Hsintien Stream.
that direction. The tidally averaged sediment concen-tration distributions, averaged over two spring-neap cycles, are shown in Fig. 13. Alocal maximum con-centration exists around the null point in the Tanshui River–Tahan Stream. The bottom sediment concentra-tions in the Tanshui River at the downriver side may be twice as much of those in the Tahan Stream at the up-river side.
6. Conclusions
Anumerical model has been developed and expended to study the hydrodynamic characteristics and suspen-ded sediment transport in the Tanshui River estuarine system. The model is a laterally averaged, two-dimen-sional, real-time model. The model has been calibrated and verified with observed time series water level, cur-rent, and salinity data and spatial distributions of tidally averaged salinity in 1994 and 1995. The residual circu-lation in the estuary has been validated by comparing model results with theoretical analysis and field data (Hsu et al., 1999). Comparisons of the suspended cohe-sive sediment concentrations were also evaluated and showed good agreement between observations and model results.
The validated model was used to investigate the ti-dally averaged salinity distributions, residual circulation and suspended sediment concentration under low flow condition. The most important aspect of estuarine cir-culation is the gravitational circir-culation. The tidally averaged circulation is characterized by the upriver movement of denser saline water in the lower layer and the downriver movement of fresher water in the upper layer. The circulation pattern is driven by the longitu-dinal gradient of salinity distribution. Under Q75 flow condition, the limits of salt intrusion are located at Hsin-Hai bridge in Tahan Stream, 26 km from the river mouth, 12 km from its mouth in the Keelung River and 4 km from its mouth in the Hsintien River. The null points are located at the heads of salt intrusion, i.e. 1 ppt isohaline, as an indicator. The tidally averaged sediment concentration distributions, averaged over two spring-neap tidal cycles, exhibits a local maximum concentra-tion around the null point in the Tanshui River–Tahan Stream.
Acknowledgements
The project under which this study is conduced is supported by National Science Council, Taiwan, under grant no. NSC 90-2211-E-002-082. The authors like to express their appreciation to the Taiwan Water Con-servancy Agency and Environmental Protection Ad-ministration for providing the prototype data.
Fig. 12. Calculated residual velocity averaged over 58 tidal cycles under Q75 flow condition in (a) Tanshui River–Tahan Stream, (b) Keelung River, and (c) Hsintien Stream.
Fig. 13. Calculated suspended sediment distributions averaged over 58 tidal cycles under Q75 flow condition in (a) Tanshui River–Tahan Stream, (b) Keelung River, and (c) Hsintien Stream.
References
Ariathurai, R., Krone, R.B., 1976. Finite element model for cohesive sediment transport. Journal of Hydraulic Engineering Division, ASCE 102, 323–338.
Barros, A.P., 1996. An evaluation of model parameterizations of sediment pathways: a case study for the Tejo estuary. Continental Shelf Research 16, 1725–1949.
De Vries, M., Klaassen, G.J., Striksma, N., 1989. On the use of movable bed models for river problems: state of art. Symposium of River Sedimentation, Beijing, China.
Delo, E.A., 1988. Estuarine muds manual. Report No. SR 164, Hydraulic Research, Wallingford, UK, 64 pp.
Dry, K.R., Evans, E.M., 1989. Dynamics of turbidity maximum in a homogeneous tidal channel. Journal of Coastal Research 5, 23– 30.
Hansen, D.V., Rattray Jr., M., 1965. Gravitational circulation in straits and estuaries. Journal of Marine Research 23, 104– 121.
Hayter, E.J., Mehta, A.J., 1986. Modelling cohesive sediment trans-port in estuarine waters. Applied Mathematical Modelling 10, 294– 303.
Hsu, M.H., Kuo, A.Y., Kuo, J.T., Liu, W.C., 1996. Study of tidal characteristics, estuarine circulation and salinity distribution in Tanshui River system (I). Technical Report No. 239, Hydrotech Research Institute, National Taiwan University, Taipei, Taiwan (in Chinese).
Hsu, M.H., Kuo, A.Y., Kuo, J.T., Liu, W.C., 1997. Study of tidal characteristics, estuarine circulation and salinity distribution in Tanshui River system (II). Technical Report No. 273, Hydrotech Research Institute, National Taiwan University, Taipei, Taiwan (in Chinese).
Hsu, M.H., Kuo, A.Y., Kuo, J.T., Liu, W.C., 1999. Procedure to calibrate and verify numerical models of estuarine hydrodynamics. Journal of Hydraulic Engineering ASCE 125, 166–182.
Hydrological Year Book of Taiwan 1996. Water Resources Planning Commission (in Chinese).
Krone, R.B., 1962. Flume studies of the transport in estuarine shoaling processes. Final Report, Hydraulic Engineering Laboratory and
Sanitary Engineering Research Laboratory. University of Califor-nia, Berkeley.
Kuo, A.Y., Nichols, M., Lewis, J., 1978. Modeling sediment move-ment in the turbidity maximum of an estuary. Bulletin 111, Virginia Water Resources Research Center. Blacksburg, Virginia. Li, Z.H., Nguyen, K.D., Brun-Cottan, J.C., Martin, J.M., 1994.
Numerical simulation of the turbidity maximum transport in Gironde estuary (France). Oceanologica Acta 17, 479–500. Mulder, H.P.J., Udink, C., 1991. Modelling of cohesive sediment
transport. Acase study: the western Scheldt estuary. In: Edge, B.L. (Ed.), Proceedings of the 22nd International Conference on coastal Engineering, ASCE, pp. 3012–3023.
McNeil, J., Taylor, C., Lick, W., 1996. Measurements of erosion of undisturbed bottom sediments with depth. Journal of Hydraulic Engineering ASCE 122, 316–324.
OÕConnor, B.A., 1971. Mathematical model for sediment distribution. Proceedings of 14th IAHR Conference, Paper D23, Paris, France. Odd, N.V.M., Owen, M.W., 1972. Atwo-layer model of mud transport in the Tames estuary. Proceedings of Institute of Civil Engineers, London, pp. 195–202.
Park, K., Kuo, A.Y., 1993. A vertical two-dimensional model of estuarine hydrodynamics and water quality. Special Report in Applied Marine Science and Ocean Engineering, No. 321, Virginia Institute of Marine Science, Gloucester Point, Virginia.
Partheniades, E., 1965. Erosion and deposition of cohesive soils. Journal of Hydraulic Engineering Division, ASCE 91, 105–139. Pritchard, D.W., 1960. The movement and mixing of contaminants in
tidal estuaries. In: Pearson, E.A. (Ed.), Waste Disposal in the Marine Environment. Pergamon, New York, pp. 512–525. Ross, M.A., Mehta, A.J., 1989. On the mechanics of lutoclines and
fluid mud. Journal of Coastal Research 5, 51–62.
Smith, T.J., Kirby, R., 1989. Generation, stabilization and dissipation of layered fine sediment suspensions. Journal of Coastal Research 5, 63–73.
Sverdrup, H.V., Johnson, M.W., Fleming, R.H., 1942. The Oceans. Prentice-Hall, New York.
Van Rijn, L.C., Rossum, H.V., Termes, P., 1990. Field verification of 2D and 3D suspended-sediment models. Journal of Hydraulic Engineering, ASCE 116, 1270–1288.
