Numerical Simulation of Circulation and Salinity
Distribution in the Tanshui Estuary
MING-HSI HSU*, ALBERT YI-SHUONG KUO**, WEN-CHENG LIU*, AND JAN-TAI KUO***
*Department of Agricultural Engineering and Hydrotech Research Institute National Taiwan University
Taipei, Taiwan, R.O.C.
**School of Marine Science / Virginia Institute of Marine Science The College of William and Mary
Gloucester Point, VA, U.S.A.
***Department of Civil Engineering and Hydrotech Research Institute National Taiwan University
Taipei, Taiwan, R.O.C.
(Received February 24, 1998; Accepted July 20, 1998)
ABSTRACT
Flow conditions in an estuary are primarily influenced by tides, seasonally varying freshwater discharge and intrusion of salt water from the ocean. A laterally averaged two-dimensional numerical model is developed and expanded for an estuarine system with tributaries as well as a mainstream. A series of comprehensive field data is used for model calibration and verification. The model simulates the time varying flow and salinity distribution in the Tanshui estuary. Reasonable agreement is obtained between the model results and observed data. This paper focuses on model experiments and applications. The model is applied to investigate the estuarine circulation, residual velocity profiles and salinity distributions under various hydrological conditions in the Tanshui estuary. The numerical experiment results show that, under Q90 flow conditions, if the wastewater is totally cut-off by wastewater interception systems, the salinity will increase by 2 to 4 parts per thousand in the middle portion of the mainstream and lower portion of the tributaries. However, the limit of salt intrusion will change little since it is limited by the shallow depth in the upper estuary.
Key Words: numerical model, residual circulation, salinity intrusion, Tanshui River system
well-mixed or vertically homogeneous (where tidal currents dominate freshwater flow) to express the relative salinity stratification (Stommel and Farmer, 1952; Cameron and Pritchard, 1963; Bowden, 1967). However, complex flow and transport processes in estuaries resulting from the interaction of tidal forcing, surface wind stress, irregular topography, and density stratification due to mixing of fresh and salt water tend to increase the residence time of pollutants in the estuary. These factors may lead to local accumulation, where pollutant concentrations reach levels that are harmful to the aquatic environment.
Numerical models designed to simulate the rel-evant transport process in estuaries should account in detail for both advective and turbulent transport. Most of the earlier numerical models of estuaries have been used to investigate tidal dynamics using the one-dimen-sional flow theorem (Harleman, 1971). The one-di-mensional dispersion model has also been used to investigate and predict salt intrusion (Thatcher and Harleman, 1972). However, application of the
one-I. Introduction
Estuaries, in the traditional sense, are regions of transition from river to ocean. They are characterized by the possibility of tidal motions communicated from the sea, and by gradients of salinity and density asso-ciated with the progressive admixture of river water and seawater. The action of gravity upon the density difference between seawater and freshwater tends to cause vertical salinity stratification and a convection flow that has come to be known as estuarine circulation. Freshwater flow and tides are the dominant variables determining flow, the distributions of salinity and circulation within the estuary. Estuaries traditionally have been classified according to their geomorphology and their salinity stratification. The terms commonly applied are : coastal plain and fjord to express the geomorphology; and salt wedge or highly stratified (where freshwater flow dominates tidal currents), par-tially mixed or moderately stratified (where freshwater flow and tidal currents are relatively balanced), and
dimensional model to a partially mixed estuary has severe restrictions. The tidal current has large vertical amplitude and phase variations, which can not be de-scribed adequately by a one-dimensional model. Salt intrusion is mainly induced by internal gravitational circulation (Pritchard, 1956), which can not be properly modeled in a vertically integrated system. The ver-tically integrated model does not appear to be very useful for circulation and water quality study in a partially mixed estuary. Often, estuaries are elongated in form so that lateral variation relatively insignificant. Dyer (1973) stated that when the ratio of river width to length is small, lateral shear may be sufficiently intense to create laterally homogeneous conditions. The estuarine reach is more than 20 km in length and the average river width is about 600 m in the Tanshui River. A laterally integrated model was developed and applied in this study to the Tanshui River estuary.
Hamilton (1975) developed a vertical two-dimen-sional numerical model of a rectangular geometry to study circulation in the Rotterdam Waterway and vertical mixing within a tidal cycle (Bowden and Hamilton, 1975). Blumberg (1975, 1977, 1978), Elliott (1976) and Rao (1995) included nonuniform geometry in their two-dimensional models of the Chesapeake Bay, Potomac River and Godavari estuary, respectively. Festa and Hansen (1976, 1978), on the other hand, applied the steady-state, uniform geometry model to study internal circulation and turbidity. Wang and Kravitz (1980) developed a semi-implicit, two-dimensional model for circulation in a partially mixed estuary. Application of their model to the Potomac River in-dicated large longitudinal and vertical changes in tide, density-driven and wind-driven circulation. Despite the fact that most of these studies were rather restricted in application, the two-dimensional model does appear to reveal several major features characteristic of a partially mixed estuary.
In this study, a laterally integrated, two-dimen-sional, real-time model of hydrodynamics and salinity is developed and expanded to handle the tributaries as well as the mainstream of an estuarine system, and is applied to the Tanshui River estuary. The Tanshui River estuary is the largest estuarine system in Taiwan, with a drainage basin that includes the city of Taipei. It consists of there major tributaries : the Tahan Stream, Hsintien Stream and Keelung River (Fig. 1). The downstream reaches of all three tributaries are influ-enced by tides. The upriver reaches are affected by daily varying freshwater discharges. The river system has a total drainage area of 2726 km2
and a total channel length of 327.6 km.
The dynamic processes involving interaction be-tween river discharge and tidal currents are complex
Fig. 1. Map of the Tanshui River system and the model transects.
and lead to a series of distinctive types of estuarine circulation in the Tanshui estuary. The numerical model is calibrated and verified using observational data. Then the model is applied to study residual circulation and salinity distribution. Simulations of the combined freshwater and tidal flow in the Tanshui estuary under various hydrological conditions have been carried out and results are reported. The model also provide further insight into residual circulation and transport processes.
II. Numerical Simulation Model
1. Governing EquationThe laterally integrated two-dimensional hydro-dynamic and transport equations are based on the principles of conservation of volume, momentum and mass. With a right-handed Cartesian coordinate system in which the x-axis is directed seaward and the z-axis directed upward, the governing equations include the laterally integrated continuity equation:
∂(uB)
∂x + ∂(wB)
∂z = qp; (1)
the cross-sectionally integrated continuity equation: ∂ ∂t(Bηη) + ∂ ∂x – H(uB)dz η = q ; (2)
the laterally integrated momentum balance equation:
∂(uB) ∂t + ∂(uBu) ∂x + ∂(uBw) ∂z
= – Bρ∂∂px +∂∂x(AxB∂∂ux) +∂∂z(AzB∂∂uz) ; (3)
the longitudinal pressure gradient results from the hydrostatic pressure assumption:
∂p ∂x=ρg ∂η ∂x + g ∂ρ ∂xdz z η ; (4)
the laterally integrated mass balance equation for salt: ∂(sB) ∂t + ∂(sBu) ∂x + ∂(sBw) ∂z =∂∂x(KxB∂∂sx) +∂∂z(KzB∂∂sz) + S0; (5)
the equation of state:
ρ=ρ0(1+ks), (6)
where t=time; η=the position of the free surface above the mean sea level; u & w=the laterally averaged ve-locities in the x- and z- directions, respectively; s=the laterally averaged salinity; B & Bη=the river width and width at the free surface, including the side storage area; H=the water depth below the mean sea level; qp=the lateral inflow per unit lateral area; q=the lateral
inflow including exchange with the side storage area; p & ρ=the pressure and water density; ρ0=the density
of freshwater; k=a constant (7.5×10–4
ppt–1
); g= gravitational acceleration; Ax & Az=turbulent
viscosi-ties in the x- and z- directions, respectively; Kx &
Kz=turbulent diffusivities in the x- and z- directions,
respectively; S0=the source or sink of salt.
2. Turbulent Closure Model
A. Vertical Turbulent Mixing Coefficients
The vertical coefficients (Az and Kz) are strongly
affected by the flow velocity, the relative roughness of the flow channel and the vertical stratification. Their values can vary over several orders of magnitude at a fixed point in an estuary during a tidal cycle (Odd and Rodger, 1978).
In the present model, the formulations for Az and
Kz are Az=αZ 2 (1 – Z h) 2 ∂u ∂z (1 +βRi) – 1/2 +αwHw 2 T exp ( – 2LπZ) (7) Kz=αZ 2 (1 – Z h) 2 ∂u ∂z (1 +βRi) – 3/2 +αwHw 2 T exp ( – 2 π L Z) , (8)
where Z=the distance from the surface; α, β and αw=constants to be determined through model
calibra-tion; Hw, T and L=the height, period and length,
respec-tively, of wind-induced waves.
Laboratory experiments indicate (Komori et al.,1983) that in a stably stratified open channel flow, vertical turbulent fluxes are closely related to the local Richardson number, defined as
Ri= – g ρ(∂ρ∂z)(∂∂uz) – 2 . (9)
B. Horizontal Turbulent Mixing Coefficients
The horizontal mixing coefficients (Ax and Kx) are
on the order of 105
cm2
/sec of the vertical mixing coefficients (Dyer, 1973). Results of diffusion mea-surements in English estuarine waters showed that Kx
ranged from 104
to 106
cm2
/sec (Talbot and Talbot, 1974). Festa and Hansen (1976) studied the important effect of exact values of Ax and Kx. Varying the
mo-mentum exchange coefficient from Ax=Az to Ax=106 Az
had negligible effects on the results of their tidal average model. The change, however, in the mass exchange coefficient from Kx=Kz to Kx=107 Kz did produce
sig-nificant change in their results.
The horizontal mixing terms, despite their relative insignificance in the momentum balance, are retained in the model for the purpose of stability. The present model uses constant values for Ax and Kx, and they are
adjusted, within the range of 104
to 106
cm2
/sec, through model calibration and verification.
3. Treatment of Branching Estuary
The model is expanded to treat the interaction between the tributaries and mainstream of an estuarine system for application to the Tanshui River system. The flow conditions at the tributary junction are solved by using expanded continuity and momentum equa-tions. The surface elevation and vertical velocity are calculated using Eqs. (1) and (2) at the mainstream-tributary confluence, with the added effect of mainstream-tributary inflows or outflows. Because of the spatially staggered grid used in the model, no representative velocity point is situated at the junction segment of the model. The momentum balance at the transects surrounding the junction segment is handled by neglecting the horizon-tal advective term ∂(uBu)∂
x and turbulent diffusion term
∂
∂x(AxB∂∂ux). The magnitudes of these two terms in the
momentum equation are assumed to be locally negli-gible in comparison with other terms (Dyer, 1973). The
mass balance equation is modified to account for the flux from or into the tributaries. The flux includes the horizontal advective term ∂(sBu)∂
x and diffusive term
∂
∂x(KxB∂∂sx). A detailed description of the treatment
between tributaries and mainstream can be found in Hsu et al. (1997).
4. Boundary Conditions
A. Free Surface
The condition of no mass flux through the free surface is effected by specifying a zero diffusion co-efficient there. The wind stress term is used to account for momentum introduced into the estuary. That is, at z=η,
ρAz∂∂uz = CDρaUw Uw , (10)
where CD is the dimensionless drag coefficient (1.3×
10–3
), ρa is the air density (1.2×10–3 g•cm–3) and Uw
is the wind speed at a height of 10 m above the sur-face.
B. Bottom
The condition of no mass flux through the bottom is effected by zero vertical velocity and diffusion co-efficient there. The bottom stress, which accounts for resistance friction at the estuarine bottom, is calculated using a quadratic law. That is, at z=−H,
τb= Az∂∂uz = Cd•u1u1 , (11)
where τb is the bed shear stress and u1 is the velocity
near the bottom. In this study, the velocity at the first grid point above the bottom is used, and Cd is the
friction coefficient on the order of 10–3
, which is:
Cd=gn2∆z–1/3, (12)
where n is the Manning’s friction coefficient and ∆z is the bottom layer thickness.
C. Upstream Boundary
The landward boundary of the model is chosen at a location upriver of the landward limit of the tidal influence. Dirichlet boundary conditions are imposed on the upstream boundary. It is assumed that the freshwater discharge, cross-sectional area and salinity are given at this boundary. The velocity at the upstream end is written as:
u =Q(t)
A , (13)
where Q(t) is the freshwater discharge through the upstream boundary and A is the cross-sectional area there, which is related to Q through the Manning for-mula. Since the tides propagates further upriver than the salt does, the salinity at the upstream boundary is specified as zero.
D. Downstream Boundary
The seaward boundary is located near the mouth of the estuary. The surface elevation is specified as a function of time with either harmonic functions or field measurements at this boundary. Since the flow condition is subcritical flow downstream, there is no need to specify the velocity at this boundary. During flood tide (u<0), oceanic water is advected into the estuary, thus increasing salinity at the mouth until oceanic salinity is reached. During ebb tide (u>0), the horizontal salinity profile is assumed to have advected out of the mouth as a “frozen” pattern, i.e., neglecting the diffusion. That is,
∂s ∂t = –∂
s
∂xu . (14)
5. Treatment of Shallow Shoals
Kuo and Park (1995) proposed a framework for coupling the shallow areas with the main channel in numerical modeling of coastal plain estuaries. They demonstrated that accounting for the mass and momen-tum exchanges between the main channel and shallow areas is essential not only for computation of the conditions in the shallow areas, but also for proper simulation of tidal propagation along the main channel. If water and momentum exchanges between the shallow areas and the main channel are not accounted for, the model can not reproduce the along-channel variations of both the tidal range and tidal phase using a single set of calibrated friction coefficients. The coupling framework proposed by Kuo and Park (1995) is refined and adopted to fit the model for the Tanshui River system. The shallow areas are treated as temporary side storage areas, thus functioning as a sink or source of momentum and mass to the main channel as the tide rises or falls, respectively.
III. Model Calibration and
Verifica-tion
A large set of accurate data is required to calibrate the numerical simulation model and to verify its
ca-Fig. 3. Daily measured river discharges at the upstream boundaries. pability to predict the flow and salinity distributions.
Sufficient data for this investigation were available from a comprehensive field study of the Tanshui River system carried out in 1994 and 1995 by the Taiwan Provincial Government Water Resources Department. The measured data were collected and analyzed for the model validation.
The numerical scheme of the simulation model is applied with a horizontal grid spacing of ∆x=1000 m and a vertical grid spacing of ∆z=2 m for the top layer and ∆z=1 m for the layers below. The thickness of the surface layer varies with time and depends on the actual calculated water level. A time step, ∆t, is limited by the Courant-Fredrick-Levy (CFL) stability condition ∆t≤ ∆x
gh . A time step increment (∆t) of 108 seconds, which guaranteed stability, was used for all the model runs.
at the upstream boundaries, respectively. These con-ditions would also serve to investigate the model re-sponse to the interaction of tidal forcing and varying river discharge. Figure 4 shows some segments of the model results, together with measured data. It shows that the upriver station (Fig. 4(b)) has much more conspicuous response to pulse of high freshwater dis-charge than the downstream station (Fig. 4(a)). An intensive survey was conducted by the Taiwan Provin-cial Government Water Resources Department on 24 June 1994. Half-hourly measurements of velocity and hourly salinity were made continuously for 13 daylight hours. Figure 5 shows the computed surface elevation and flow discharge in one tidal cycle, together with time series measured data. The comparison shows that the model can faithfully reproduce tidal propagation, tidal flow and river flow. All available velocity data for the Tanshui River system were obtained during only 13 hours in this survey. No time varying residual velocity could be deduced through standard low-pass filtering to discern if two-layer estuarine circulation Fig. 4. A comparison of the computed surface elevation with field data. (a) Tu-Ti-Kung-Pi, (b) Hsin-Hai Bridge.
Manning’s friction coefficient is the most impor-tant calibration parameter affecting the calculation of the surface elevation, velocity and flow. Since tidal flow constitutes the major portion of energy in estua-rine flows, the Manning’s coefficient was adjusted based on comparison of the predicted tidal wave propa-gation with measured data. Both the spatially varying tidal range and phase were calibrated with results presented in Hsu et al. (1977). The calibrated model has a Manning’s friction coefficient of from 0.032 to 0.026 for the Tanshui River-Tahan Stream, 0.015 for the Hsintien Stream, and from 0.023 to 0.016 for the Keelung River.
The friction coefficient was verified through simu-lation of prototype conditions during the period March 15, 1994, to September 30, 1994. Hourly measure-ments of the water surface elevation at the river mouth and daily freshwater discharges upriver of the tidal limits were used as boundary conditions. Figures 2 and 3 show portions of the measured water surface eleva-tion at the river mouth and daily freshwater discharges
turbulent diffusion coefficients. The turbulent mixing term in the vertical direction is the dominant factor which determines stratification in the water column. The same prototype simulation for friction coefficient verification was used for calibration of constants in the Fig. 5. Model simulation and field measurement of the surface elevation and flow on 24 June 1994. (a) River Mouth, (b) Kuan-Du Bridge,
(c) Taipei Bridge, (d) Hsin-Hai Bridge, (e) Pa-Ling Bridge.
exists in the Tanshui River system. However, theo-retical analysis indicated that two-layer estuarine cir-culation may exist in the lower portion of the Tanshui estuary (Hsu et al., 1998).
Fig. 6. A comparison of time series data between the computed and measured salinity on 24 June 1994. (a) Kuan-Du Bridge (top layer, middle layer and bottom layer represent 1.0 m, 5.5 m and 10.5 m below the water surface, respectively); (b) Taipei Bridge (top layer and bottom layer represent 1.0 m and 7.5 m below the water surface, respectively); (c) Hsin-Hai Bridge (top layer and bottom layer represent 1.0 m and 2.5 m below the water surface, respectively); (d) Pa-Ling Bridge (top layer and bottom layer represent 1.0 m and 6.5 m below the water surface, respectively).
as constant for the main stem and each of the tributaries. The calibrated values for the longitudinal dispersion coefficient are Ax=Kx=28×105 cm2/s for the Tanshui
River-Tahan Stream, Ax=Kx=3.5×105 cm2/s for the
Hsintien Stream, and Ax=Kx=5×104 cm2/s for the Keelung
River. Since the wind was weak on June 24, 1994, no wind was included in the model simulation for that period. Figure 6 shows a comparison of the time series salinity distributions between the computed and mea-sured data in different layers. The results show that the numerical model can favorably determine the trend of salinity distribution. The constants in the turbulent closure model were verified by simulation of prototype conditions in 1995 (Hsu et al., 1997).
Both the field data and model simulation (Fig. 6) indicate that the vertical salinity difference varies with the tidal phase from nearly zero to as much as 6 ppt turbulent closure model. The maximum salinities
mea-sured at the river mouth on each survey date were linearly interpolated in time and used for the boundary condition, which specifies the high tide salinity during each tidal cycle.
The model has daily wind speed and direction as input data and computes the wind speed in the longi-tudinal direction at each segment. The longilongi-tudinal wind speed is then used to compute the surface wind stress. The formulation for the vertical diffusion co-efficient does not explicitly account for the wind di-rection. The wind speed and direction are accounted for through specification of period, height and length of wind wave in Eqs. (7) and (8).
The calibrated constants for vertical mixing are α=0.0115 and β=0.75. The coefficient for wind mixing (αw) was calibrated to be 5×10–2. They are specified
Fig. 8. Computed salinity distributions at spring tide under mean flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream. (The numbers on the contours refer to the salinity in parts per thousand)
Table 1. The Amplitudes and Phases of Nine Tidal Constituents at the River Mouth
Constituents Amplitude Phase
(ak, cm) (φk, degree) M2 104.87 35.60 S2 27.85 –7.18 N2 21.33 –43.76 K1 19.94 –133.79 Sa 17.06 –143.32 O1 16.22 –47.96 K2 7.17 131.50 P1 7.14 –110.58 M4 2.83 47.44
Table 2. Freshwater Discharges Imposed under Upstream Boundary Conditions and Lateral Inflow
Upstream freshwater Mean flow condition Q90 flow condition
discharges/Lateral inflow
Tahan Stream 62.1 m3/s 4.0 m3/s Hsintien Stream 72.7 m3/s 6.9 m3/s Keelung River 26.1 m3/s 1.3 m3/s Lateral inflow from
51.13 m3/s 2.89 m3/s watershed
Lateral inflow from
20.63 m3/s 20.63 m3/s a wastewater
a
This flow is assumed to be zero with wastewater interception.
Fig. 7. The time varying water surface elevation at the river mouth, showing two spring-neap tidal cycles in September.
at Kuan Du. Hsu et al. (1997) also showed that the spatial distribution of the tidal average salinity has a vertical difference of one to two parts per thousand. Therefore, the Tanshui estuary may be classified as a partially mixed estuary.
IV. Numerical Model Experiments
The calibrated and verified model was used to perform a simulation of combined freshwater and tidal flows under various hydrological conditions to inves-tigate the flow and salinity structure as well as the estuarine circulation in the Tanshui River system. The values of all the calibrated constants were determined as described in a previous section, and no further ad-justment was made. Nine-constituent tide, extracted from harmonic analysis, was used to specify the down-stream boundary conditions for all the model runs. We performed Fourier analysis of the surface elevation ηFig. 9. Computed salinity distributions at neap tide under mean flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
where ak and φk are the amplitude and phase angle of
the k-th harmonic, respectively; η is the mean value of the tidal elevation; σk is the angular velocity. The
historical long-term mean value of the tidal elevation is 6.5 cm. Table 1 lists the amplitudes and phases of nine tidal constituents at the river mouth, which were derived from harmonic analysis of the field data.
Three hydrological conditions were investigated: mean freshwater discharge, Q90 discharge and Q90 with
wastewater interception (Table 2). The table also shows the lateral inflow from the watershed below the gauging stations and from wastewater. Historical long-term mean discharges were used for the upriver boundary conditions in the tributaries in the first simulation. The mean discharges at the tidal limits of the three major tributaries are 62.1 m3
/s, 72.7 m3
/s and 26.1 m3
/s for the Tahan Stream, Hsintien Stream and Keelung River, respectively. A high tide salinity of 25 ppt at the river mouth was used for the mean flow simulation. The
model simulated one-year (705 tidal cycles) hydraulic and salinity conditions. Figure 7 shows an example of the water surface elevation for two spring-neap tidal cycles at the river mouth. The average velocity and
Fig. 10. Computed residual velocity profiles at spring tide under mean flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
tide provides more tidal mixing and dispersion, leading to upstream movement of saline water. The most stratification occurs 7 to 10 km from the river mouth in the Tanshui River-Tahan Stream. There is little salt intrusion in the Hsintien Stream. Figure 10 shows the residual current at spring tide. Estuarine circulation occurs only in deeper waters around Kuan-Du under mean flow conditions. An important capability of the coupled hydrodynamic and salinity model is its ability to calculate the residual circulation. The residual cir-culation arises from the pressure gradient generated by the longitudinal salinity variation, a baroclinic force which increases with the depth as well as the salinity gradient. Consequently, the residual circulation tends to be strongest in deep sections of the estuary and in regions where the salinity gradient is largest. This expectation is borne out by the numerical calculations Fig. 11. Schematic diagram of the residual circulation.
Fig. 12. Computed salinity distributions at spring tide under Q90 flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream. (The numbers on the contours refer to the salinity in parts per thousand)
Fig. 13. Computed salinity distributions at neap tide under Q90 flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
salinity over a three day period were used to represent the residual current and tidal average salinity during spring or neap tide. Figures 8 and 9 illustrate the tidal average salinity during spring and neap tides, respec-tively. The limit of salt intrusion is represented by 1 ppt isohaline. It shows that the spring tide pushes the limit of salt intrusion farther upriver in the Tanshui River-Tahan Stream and Keelung River. The spring
Fig. 14. Computed residual velocity profiles at spring tide under Q90 flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
as shown in Figs. 8(a) and 9(a). The present results indicate that the large variations in the depth charac-teristic of the downriver section of the Tanshui estuary cause recirculating regions to develop as shown in
Fig. 11.
To study the responses of residual velocity and salinity distribution to varying river discharge, another model simulation was conducted using Q90 flow. (Q90
is a flow with an exceedance probability of 90%.) Q90
represents a very low flow condition but has an occurrance probability of 10%. The scatter plots of field data reveal no definite relation between salinity in the river mouth and freshwater discharge upriver. Thirty-two ppt of high tide salinity at the river mouth, obtained and estimated from the available field data, were used in the simulation. All other conditions were kept the same as in the mean flow simulation. Figures 12 and 13 present salinity distributions during spring and neap tides, respectively. They show that Q90 flow
conditions push the limit of salt intrusion farther upriver than do mean flow conditions. The limit of salt in-trusion is located at the Hsin-Hai bridge in the Tahan Stream, 14 km from the Keelung River mouth in the Keelung River, and at the Chung-Cheng bridge in the Hsintien Stream. This extensive intrusion of saline water introduces a significant baroclinic force and induces a strong residual circulatory system in the estuary. Figures 14 and 15 present residual velocity profiles during spring and neap tides, respectively. They indicate that estuarine circulation exists in all three tributaries as well as in the mainstream Tanshui River. Figures 14 and 15 also show that the spring tide provides more turbulent mixing energy and, thus, has weaker residual circulation than does the neap tide.
Comparison of two model simulations described above does not preclude the possibility that the increase in the salt intrusion under Q90 conditions may result
from an increase in salinity boundary conditions at the river mouth. Therefore, a model sensitivity run was conducted using 32 ppt at the river mouth in the mean flow simulation. The results (not shown) show little movement of salt intrusion limits based on the mean flow simulation with 25 ppt salinity as boundary con-ditions. The increase in salinity at the downstream boundary has little effect in the upper reaches of the estuary. Therefore, the freshwater flow reduction is the dominant factor pushing saline water upriver.
The third model simulation was run for one-year (705 tidal cycles) with Q90 flow as the upstream
bound-ary condition while the lateral inflow from wastewater along the river was totally cut-off by the implemented interception system. Due to the diversion of a large amount of wastewater to the ocean outfall system, the freshwater flow in the river significantly decreased. Figures 16 and 17 show salinity distributions at spring and neap tides, respectively. The limit of salt intrusion is located near the Hsin-Hai bridge in the Tahan Stream,
Fig. 15. Computed residual velocity profiles at neap tide under Q90 flow conditions. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
at a 15 km distance from the Keelung River mouth in the Keelung River, and at a 7 km distance from the Hsintien Stream mouth in the Hsintien Stream. Com-parison of the simulation results between no
wastewa-Fig. 16. Computed salinity distributions at spring tide under Q90 flow conditions while wastewater was totally cut-off. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream. (The numbers on the contours refer to the salinity in parts per thousand)
ter cut-off (shown in Figs. 12 and 13) and total wastewater cut-off indicates that salinity increases about 2 ppt in most of the estuary as a result of wastewater diversion. The greatest increase occurs in the Keelung River with an increase of 4 ppt around 5 km from the Keelung River mouth. However, the limit of salt in-trusion moves only slightly upriver if the wastewater is totally cut-off. The model results indicate that the limits of salt intrusion stay near topographic sills. The sills at the upriver end play an important role in pre-venting further salt intrusion as freshwater input dimin-ishes.
V. Conclusions
Using a vertical (laterally averaged) two-dimen-sional numerical model and calibrated values for model
Fig. 17. Computed salinity distributions at neap tide under Q90 flow conditions while the wastewater was totally cut-off. (a) Tanshui River-Tahan Stream, (b) Keelung River, (c) Hsintien Stream.
constants, analysis has been carried out on the com-puted residual circulation and salinity distribution in the Tanshui estuary under both mean and Q90
fresh-water discharge conditions. The following conclusions have been made:
(1) Both the field data and model simulation results indicate that the Tanshui estuary system may be classified as a partially mixed estuary. (2) The current velocity in estuaries may be
decom-posed into two components, tidal and residual. In a partially-mixed estuary, the dominant re-sidual velocity is characterized by upriver move-ment of more saline water in the lower layer and downriver movement of fresher water in the upper layer. The computed residual velocity profiles in the Tanshui estuary show these two-layer flow characteristics, particularly under low river discharge conditions.
(3) With given river freshwater discharge condi-tions, the computed salinity distributions and residual velocity profiles vary over the spring-neap tidal cycle. The results show that the spring tide leads to more turbulent mixing and, thus, has less vertical stratification and weaker re-sidual circulation. The enhanced mixing during spring tide also increases the salt intrusion. (4) The magnitude of river discharge is the dominant
factor affecting salinity distribution in the Tanshui estuary system. Under mean flow conditions, sea water barely intrudes into the Hsintien Stream. In the Keelung River, the limit of salt intrusion varies from 9 to 5 km over the spring-neap cycle. Reduction of river discharge to Q90 flow increases the extent of salt intrusion
by about 6 km, 6 km and 7 km in the Tahan Stream, Keelung River and Hsintien Stream, respectively.
(5) The total diversion of wastewater under Q90 flow
conditions results in only a slight change in the limit of salt intrusion since it is limited by the shallow depth in the upriver regions. However, there is significant salinity increase in the middle reach of the Tanshui River and in the lower portions of all the tributaries. The most severe impact is noted in the Keelung River, where a salinity increase as high as 4 ppt is predicted. This may be attributed to the more gentle bottom slope in this tributary.
Acknowledgment
The project under which this study was conducted, in part, was supported by the National Science Council, R.O.C., under grant No. NSC 86-2621-E-002-013. This financial support is highly ap-preciated. The authors thank the Taiwan Provincial Government Water Resources Department for providing prototype data. The authors also would like to express their appreciation to the manu-script reviewers, whose comments led to substantial improvement of this paper.
References
Blumberg, A. F. (1975) A Numerical Investigation into the Dynamics
of Estuarine Circulation, p. 59. Chesapeake Bay Institute Tech.
Rep. 91, The Johns Hopkins University, Baltimore, MD, U.S.A. Blumberg, A. F. (1977) Numerical model of estuarine circulation.
Journal of Hydraulic Division, ASCE, 103, 295-319.
Blumberg, A. F. (1978) The influence of density variations on estuarine tides and circulation. Estuarine and Coastal Marine
Science, 6, 209-215.
Bowden, K. F. (1967) Circulation and diffusion. In: Estuaries, pp.15-36. G.H. Lauff Ed. American Associate for the Advance-ment of Science, Publication No.83, AAAS, Washington, D.C., U.S.A.
Bowden, K. F. and P. Hamilton (1975) Some experiments with a numerical model of circulation and mixing in a tidal estuary.
Estuarine and Coastal Marine Science, 3, 281-301.
Cameron, W. M. and D. W. Pritchard (1963) Estuaries, pp. 306-324. In: The Sea, Vol. 2. M.N. Hill Ed. Interscience, New York, NY, U.S.A.
Dyer, K. R. (1973) Estuaries: a Physical Introduction, p. 140. J. Wiley & Sons, New York, NY, U.S.A.
Elliott, A. J. (1976) A Numerical Model of the Internal Circulation
in a Branching Estuary, p. 85. Chesapeake Bay Institute, Spec.
Rep. 54, The Johns Hopkins University, Baltimore, MD, U.S.A. Festa, J. F. and D. V. Hansen (1976) A two-dimensional numerical model of estuarine circulation: the effects of altering depth and river discharge. Estuarine and Coastal Marine Science, 4, 309-323.
Festa, J. F. and D. V. Hansen (1978) Turbidity maxima in partially mixed estuaries: a two-dimensional numerical model. Estuarine
and Coastal Marine Science, 7, 347-359.
Hamilton, P. (1975) A numerical model of vertical circulation of tidal estuaries and its application to the Rotterdam waterway.
Geophys. J. Roy Astron. Soc., 40, 1-12.
Harleman, D. R. F. (1971) One-Dimensional Model. Estuarine
Modeling: An Assessment, pp. 34-85. TRACOR, Inc., New York,
NY, U.S.A.
Hsu, M. H., A. Y. Kuo, J. T. Kuo, and W. C. Liu (1997) Study of
Tidal Characteristics, Estuarine Circulation and Salinity Dis-tribution in Tanshui River System. Technical Report No. 273,
Hydraulic Research Lab., National Taiwan University, Taipei, Taiwan, R.O.C.
Hsu, M. H., A. Y. Kuo, J. T. Kuo, and W. C. Liu (1998) Procedure to calibrate and verify numerical models of estuarine
hydrody-namics. J. of Hydraulic Engineering, ASCE, (in press). Komori, S., H. Ueda, F. Ogino, and T. Mizushina (1983) Turbulence
structure in stably stratified open-channel flow. Journal Fluid
Mechanics, 130, 13-26.
Kuo, A. Y. and K. Park (1995) A framework for coupling shoals and shallow embayments with main channels in numerical mod-eling of coastal plain estuaries. Estuaries, 18(2), 341-350. Odd, N. V. M. and J. G. Rodger (1978) Vertical mixing in stratified
tidal flows. J. of the Hydraulic Division, ASCE, 104(HY3), 337-351.
Pritchard, D. W. (1956) The dynamic structure of a coastal plain estuary. Journal of Marine Research, 15, 33-42.
Rao, A. D. (1995) A numerical modelling study of the flow and salinity structure in the Godavari estuary, east coast of India.
International Journal for Numerical Methods in Fluids, 21,
35-48.
Stommel, H. and H. G. Farmer (1952) Abrupt change in width in two-layer open channel flow. Journal of Marine Research, 11, 205-214.
Talbot, J. W. and G. A. Talbot (1974) Diffusion in shallow seas and in English coastal and estuarine waters. Rapp. P. Reun.
Cond. Int. Explor. Mer., 167, 93-110.
Thatcher, M. L. and D. R. F. Harleman (1972) A Mathematical Model
for Prediction of Unsteady Salinity Intrusion in Estuaries, p. 232.
Report No. 72-7, MIT, Cambridge, MA, U.S.A.
Wang, D. P. and D. W. Kravitz (1980) A semi-implicit two-dimen-sional model of estuarine circulation. Journal of Physical
Ocean-ography, 10(3), 441-454.
