+Bioline International Official Site (site up-dated regularly)

International Journal of Environmental Research, Vol. 3, No. 4, 2009, pp. 671-680

School of Civil Eng., Iran University of Science and Technology, P.O. Box 16765163, Tehran, Iran
*Corresponding author E-mail: etemad@iust.ac.ir

The Karkheh Dam Reservoir with a capacity of more than 5 billion cubic meter is the largest dam in Iran with both agricultural and drinking usages. Its hydrodynamics and water quality were modeled and simulated to analyze the total maximum daily load (TMDL) of Total Dissolved Solids (TDS). The simulation was supported with measurements of temperature and TDS measurements during two years. A laterally averaged 2D model called CE-QUAL-W2 was used for the simulation and hypothetical low height spillways were implemented in the model to avoid drying of the cells in the river branch. The model was then calibrated successfully with an absolute mean error of 0.71^° C. More importantly, vertical stratification of temperature and TDS in the Karkheh Reservoir was reproduced by the model throughout years 2000 to 2003. The calibrated model was then used to simulate water quality response to various TDS reduction scenarios. Model results reveal that a 50% reduction of the TDS load is required for a 40% reduction of TDS in the reservoir outlet. The modeling of a complex combination of a steep and long river –reservoir system was another important achievement of this study.

Key words: Total dissolved solids, Simulation, CE-QUAL-W2, Water quality, Modeling, Stratification

Damming of rivers has a major impact on the natural water resources and impoundments change the characteristics of a water body. This change is not limited to the hydrology but also affects the physical, chemical and even biological characteristics of water bodies (Friedl and Wuest 2002). The changes are most prominently local and water quality models can be useful in evaluation of these effects before, during and even after construction of dams. In addition, water authorities have required a Total Maximum Daily Load (TMDL) analyses to determine the necessary water management actions for elimination of the water quality impairments. Nowadays, models of watersheds and receiving water bodies are often an integral part of the total maximum daily load TMDL process. The development of a TMDL is necessary for important water bodies. The TMDL process requires the determination of the point source load and nonpoint source load allocations for a water body that is necessary to meet specified water quality objectives (De Pinto et al., 2004). Water quality models may be considered as one of the best tools available for determining the quantitative relationship between loads and water body response. Furthermore, models can be used to estimate watershed loads for existing conditions, and evaluate the effectiveness of proposed control alternatives in reducing loads and improving water quality to meet standards (Bowen and Hieronymus 2003).

The Karkheh Dam Reservoir with a capacity of more than 5 billion cubic meter is the largest dam in Iran. The importance of this multipurpose (agricultural and drinking) water resource and its high TDS urge us to model temperature and salinity in this river-reservoir system. This paper reports on the calibration, verification and application of a two-dimensional laterally averaged model for Karkheh Dam Reservoir to support a TMDL analysis of the watershed. The objective of the TMDL study is to predict the required watershed Total Dissolve Solids (TDS) load reduction to meet the water quality standard in the reservoir outlet. This paper is organized as follows: The next section will give a brief description of the study area. Then, the model formulation and application is discussed. Next, the results obtained from the calibration and verification procedures are presented. This is followed by the application of the model for TMDL analysis. The summary and conclusions are presented in the last section.

Watershed area of Karkheh River is more than 50,000 square kilometer in the west–south of Iran (fig. 1) and the dam on the River is located at 32º 27´ N and 48º 8´.

Due to the existence of Zagros Mountains at the north and east of the river, Karkheh water shed has a varying climate. The annual average rainfall varies from 300 to 800 mms, more than 50% of which occurs in winter and the rest mainly in spring and fall. The upper watershed has an area of 42644 square km with an annual average discharge of 187.6 m³s^-1, annual mean temperature of 24.6^°C and an annual mean humidity of 45.5%. The core of the dam is made of clay and its crest height is 124 meters. The dam^’s crest length is 3035 meters and is 234 ms above the sea level. Demanded water from Karkheh reservoir is supplied by the Power station outlet, water supply gates and the Dashte abbas tunnel.

CE-QUAL-W2 (Cole and Wells 2001) is a twodimensional laterally averaged hydrodynamic and water quality model capable of modeling stratified water bodies with interconnected rivers, reservoirs and estuaries (Garvey et al., 1998; Gunduz et al., 1998; Kurup et al., 1998; Lung and Bai, 2003; Kuo et al., 2006; Liu et al., 2006; Zahed et al., 2008). A major feature of the model is its ability to calculate the two-dimensional velocity field for narrow stratified water bodies. This code is developed by the waterways experiment station of corps of engineers and allows the user to include hydraulic structures such as pipes, spillways and gates in the system. Its water quality module can simulate 21 constituents.

In the model applications, the hydrodynamic runs provide real-time simulations of velocities, temperature, and a conservative tracer such as salinity prior to the water quality calculations. The model uses a numerical scheme for a direct coupling between hydrodynamic and water quality simulations. CE-QUAL-W2 also uses the same spatial grid for hydrodynamic and water quality. CE-QUAL-W2 is based on a finite-difference approximation to the laterally averaged equations of fluid motion including: the free surface wave equation; hydrostatic pressure; horizontal momentum; continuity; constituent transport and equation of state. The model quantifies the free surface elevation, pressure, density, horizontal and vertical velocities, and constituent concentrations.

In a Cartesian coordinate system with x axis directed seaward and z axis directed downward. The governing equations for an incompressible flow are as follows: The equation of continuity: where, U and W are the laterally averaged horizontal and vertical velocities, respectively, B is the width of the layer and q is the lateral inflow per unit width (Cole and Wells, 2001).

The free surface equation: where η is the free surface position, t is time and h is the total depth.

The x and z momentum balance equations, using the hydrostatic approximation, are:

where g is the gravitational acceleration, a is the channel angle, t_xx and t_zz are the turbulent shear stresses on x and z faces of the cell acting on x direction, respectively. P is the pressure and ρ is the density. The equation of state in the model is:

where T is the water temperature, TDS is the concentration of total dissolved solids and SS is the concentration of suspended solids. Finally, the equation forconservation of constituent concentration (or heat) is (Cole and Wells, 2001):

where s is the laterally averaged concentration , D_x and D_z are the horizontal and vertical eddy diffusivity, respectively and S is the lateral source/ sink term.

By solving the free surface elevation implicitly, the restriction of the Courant surface gravity wave stability criterion is lifted. Therefore, longer time steps can be used for efficient computations. Explicit numerical schemes are also used to compute velocities, which affect the transport of energy and biological/chemical constituents.The time increment is calculated automatically in the model using the following criteria (Cole and Wells, 2001).

where, ∆t is the time step, Ax is the horizontal eddy diffusivity, Az is the vertical eddy diffusivity, ∆x is the cell length, ∆z is the cell thickness, Q is the cell total discharge, V is the cell volume and ∆ρ is the density difference. In each time step, first the water levels are calculated in numerical solution. Having the new water levels, horizontal and vertical velocities are calculated and then the new 1000 concentrations are computed. Using the new horizontal and vertical velocities, the water level is calculated simultaneously (Cole and Wells, 2001).

The version 3.2 of the model provides a numerical scheme, ULTIMATE, in the advection term of the mass transport equation to eliminate the numerical dispersion and oscillation (Cole and Wells, 2001). This is very important to reproduce the vertical salinity stratification in the water column by accurately quantifying the vertical advection mass transport.To estimate the rate of the vertical eddy diffusivity/viscosity, five formulations can be used in the model. These formulations are listed in Table 1. In this Table, l_mis the mixing length, Ri is the Richardson number, k is the Von-Karman constant, u_* is the shear velocity, C is constant (0.15), t_wyis the wind shear, n ?is the molecular viscosity and C₁ is an empirical constant (100). Y(x)= max (x,0).

MATERIALS & METHODS

The structures which are modeled are: an overflow, two outlets (the bottom outlet and the power station outlet) and Dashte abbas tunnel that is considered as a selective withdrawal. The modeling period is from Feb., 19, 2001 to Dec., 1, 2003. The start time is at 24:00 o’clock on Jan, 1, 2001.Meteorological and hydrological data were obtained from the daily reports of Karkheh Meteorology Station. Hourly data of meteorology and daily data of hydrology are used in the river– reservoir modeling. Karkheh River^’s discharge, measured from 2001 to 2003 is shown in fig. 2. As seen, from end of May to early November only the base flow exists and in other periods, i.e. spring, the surface flow due to rainfall peaks up to 900 cms. Reservoir^’s measured water levels during 2000-2003 are also shown in fig. 3. The general trend is a periodic increase in the storage volume with some local decrease during dry months. These data were used for calibration of the model.

Fig.4 shows 104 kilometers of Karkheh riverreservoir system extending from Jelogeer station to the dam. In order to model this river-reservoir system, it was divided into two individual water bodies and three branches.The first water body involves two branches: the river branch and river

The second water body includes a branch of the reservoir with a zero bed slope. This water body is divided into 62 layers. The thicknesses of the layers in this water body are varying as follows: The bottom layers thicknesses vary from 5 meters at level 113 to 2 meters thick at level 150 MSL. Thermocline may exist in the layers above level 150 meters. Therefore, the layers thickness in the thermocline proximity were selected to be 1.5 meters thick (Table 2).

As shown in Table 2 and fig.4. the system is divided into river branch, river-reservoir branch and reservoir branch. The first branch is a part of river which is 220 meters above the sea level. This branch starts from kilometer zero (Jelogeer station) and extends to kilometer 49 (level 220). The second branch is a part of river-reservoir between 160 and 220 meters levels. It starts from kilometer 49 of the river (level 220) and extends to kilometer 77 (level 160). This branch is divided into 38 segments of 1 km length and a mean slope of 0.15%. It is separated from the reservoir due to the large variation of water level in the reservoir, which makes it behave like a reservoir and/or river depending on the water surface level in the reservoir. In fact, the lowest part of the reservoir (level 160) is the separation point between river and reservoir.

The third branch is the reservoir from level 113 m to level 160 m. This branch starts from km 77 of the reservoir (level 160 m) and extends to km 104 of the river (dam location). At the end of this branch there are two outlets: the power station outlet and the bottom outlet. There is also a tunnel called Dashte abbas at segment 88 of this branch at level of 176.5, which has been modeled as a selective withdraw structure. It should be added that the first simulations of the system was not successful. This was mainly due to the large difference between bottom levels in the upstream and downstream parts of the system that results in drying of some cells in the river branch. After several manipulations, this problem was overcome by using two low height artificial spillways at the end of the river branch and reservoir.

To calibrate the water level, the data from 2000 to 2001 were used. For calibration purposes, model parameters were adjusted to minimize the error between observed data and simulated ones. The water level in year 2001 was simulated with an average mean error (AME= ∑ ¦ observed -measured/N) of 0.51 m (see fig.5). Before simulating the temperature profile, light extinction coefficient was obtained from measured depth of Secchi disk, using the following equation (Williams et al., 1981):

λ= 1.11Z_s^-0.73(8)

Where, λ is the light extinction coefficient (m^-1) and Z_s is the depth of Secchi disk (m).The extinction coefficient of 0.3 was used in this study.

CE-QUAL W2 has five formulations for estimation of vertical eddy diffusivity (Table 1). Two of them, i.e. W2N and W2 consider wind effects and are recommended for deep and stratified water bodies. Similar to the other studies, our simulation revealed that W2N performs better than W2 and it was used for further steps. The major difference between W2N and W2 is that in W2N, in contrast to W2, the mixing length is grid size independent (see Table 1). For verification of the model, data from year 2003 were used and the model was run from 2001 to 2003. The simulated water level for this period in these years is shown in fig. 6. As is clear seen, the model mimics the measurements closely with an AME of 0.59 m.

Water temperature profiles were used to evaluate the model’s skill. Temperature structure was calibrated using year 2001-2002 data and was verified using year 2003 data (figs. 7 and 8). As seen, the model has reproduced the depth of surface layer and thermocline in different months very well. The AME of temperature simulation was 0.71^° C and the maximum error was 0.94^°C.For modeling of TDS, the model was run using available data from year 2003. A comparison between measured and modeled TDS is shown in figure 9.Here, the AME varied between 29.68 to 31.47 mg/l. Further studies showed that the TDS value at hypolimnion becomes maxima in late spring and early summer while the reservoir becomes homogenous during winter. Fig.9 implies the existence of an underflow that is mainly due to denser river water flowing into the reservoir.

The strong stratification was seen in fall with a maximum difference of TDS values between surface and bottom layers.Finally, in order to study the future variations of TDS and TMDL of the river-reservoir system to apply a suitable management measure, different scenarios were considered (constant, decreasing and increasing annual mean values of TDS concentrations in the river). In all of them, the following assumptions were made:(a) A three years averaged hydrological and meteorological data were used for forcing the system, (b) a 5 days time step was used to reduce the computational time, (c) the TDS load was varied linearly in time (if required). Fig. 10 displays the time series of TDS at different locations assuming no change in the annual mean value of TDS. The observed variation of TDS in the upstream is due to the changes in river discharge; high TDS values are associated with low discharge (dry seasons) and vice versa. The observed reduction of TDS concentration in the reservoir is mainly due to the dilution while the time lag between peaks of TDS values in the upstream and intakes are due to their large horizontal distance.

This figure also shows that TDS in the Dashte abbas intake and the power plant outlet will be less than the standard for agricultural use (1000 mg/L) but an increase of 10% in the TDS of Karkheh River will result in water quality standard violation during fall. In addition, it can be conferred that a 50% reduction of TMDL of TDS is required to comply with the Iran’s standard for TDS in drinking water. It should be mentioned that the Karkheh River has a maximum TDS of 1350 mg/L during September while the maximum TDS value at reservoir is 25 percent less than that and happens nearly one month later. One of the restrictions of this model is that when the depth of the water in any segment becomes less than the thickness of two layers, the model encounters dry cell and fails. In order to overcome this difficulty, layers with 0.5 m thicknesses were used for simulation of river branch. The model reduces the time step automatically to prevent rapid changes in water level. Therefore, the execution time increases significantly. In our case, the simulation of one year period took about two hours using a Pentium 4.

The study presents the application of water quality and hydrodynamic models for Karkheh River-Reservoir system in Iran.A twodimensional hydrodynamic and water quality model, CE-QUAL-W2, from the U.S. Army Corps of Engineers Waterways Experiment Station was set up. The model was calibrated and verified with data from 2000-2003 in the Reservoir. In reservoir modeling of the system without river, the incoming river water will normally be completely mixed with the surface layer of reservoir and produces inappropriate vertical structure. This was prevented by appropriate modeling of river-reservoir system. Furthermore, artificial low height spillways were implemented in the model to avoid drying of the cells in the river branch and crashing of the simulation. A close match was produced between the simulated and measured water surface elevation for model calibration and verification in the reservoir. Thermal stratification was also reproduced with quite reasonable prediction of the depths of the thermocline during the simulation period. The AME of temperature simulation was 0.74° C when using W2N scheme. The model also mimics spatial and temporal concentration distributions of key water quality constituents such as TDS. The comparison between calculated results and field data are reasonably consistent with an AME of 28 mg/L. The calibrated model was used to evaluate watershed management practices to manage the water quality in reservoirs. Results of the model scenario runs showed that a 50% reduction of the TDS loads will improve the water quality.It is believed that the used model provides a useful tool to predict the improvement of reservoir water quality resulting from various management strategies in this watershed.

The first author was the recipient of a University of Western Australia senior Gledden fellowship. We are grateful to Ministry of power, water resources management Organization for their support and providing the field data. We are also grateful to Dr. Tom Cole for his invaluable support and providing the model. We thank Dr. Se-Woong Chung for his fruitful comments on the earlier draft of this manuscript. This work was partly supported by the deputy of research, Iran Uni. of Science and Technology.

MATERIALS & METHODS

The following images related to this document are available:

Photo images