Application of a 3d Layer Integrated Numerical Model of Flow and Sediment Transport Processes to a Reservoir

Details are given in the paper of the refinement of a three-dimensional layer integrated turbulence model and its application to a scaled physical model of a reservoir, named Hamidieh Reservoir, in Iran. The strong turbulent flows generated in this reservoir are due to the high volumes of flow diversion, with low head differences. In this paper, a refined numerical model is applied to a reservoir, associated with a dam, water intakes and sluice gates, with the aim being to investigate the flow patterns and sediment transport processes in the vicinity of such hydraulic structures. The calibration of the model is undertaken using measurements made from a scaled physical model. The numerical model is able to predict the conditions for a range of scenarios that are difficult to test in the physical model. Different scenarios are introduced to investigate the effects of various intake and sluice gate configurations, as well as their operation schemes on the flow and sediment transport processes into Hamidieh Reservoir. The results indicate that the flow velocity field in the vicinity of one of the intakes varies significantly. Moreover, the sluice gates do not appear to have any considerable effect on the suspended sediment concentrations moving through the intakes.


Introduction
When a reservoir forms in a river section due to the construction of a dam, the natural system's hydro-morphological conditions will generally change significantly.Two major hydraulic changes generally occur with the construction of such a dam.Firstly, the water area upstream of the dam will generally change from lotic (i.e., running water) to lentic (i.e., still water), resulting in changes in the hydrologic, hydraulic, sediment transport and ecological processes.Secondly, diurnal and seasonal variations in the demand for water will cause short-term and long-term variations in the discharge rate.
A good understanding of the flow behaviour and the sediment transport process in reservoirs, especially near hydraulic structures, such as sluice gates and intakes, is necessary in order to design appropriate operating strategies to reduce sedimentation in the reservoir and irrigation canals, as considered herein.
The primary objective of this research programme is therefore to study the hydrodynamic and sediment transport processes in the vicinity of the hydraulic structures of a reservoir.Fluid flows in nature are three dimensional and usually highly turbulent.The three dimensionality of the turbulence and sediment transport processes in the vicinity of hydraulic structures generally creates complex conditions.A number of numerical models, including two-dimensional (2D) depth-integrated (e.g., Galppatti and Vreugdenhil [1]; Celic and Rodi [2]) and three-dimensional (3D) sediment transport models (e.g., O'Connor and Nicholson [3]; Lin and Falconer [4]), have been developed to simulate these transport processes.Van Rijn [5] and van Rijn et al. [6] investigated the influence of the basic physical parameters and developed a model for suspended sediment transport in gradually varying steady flows, in which the flow velocities were computed using a 2D depth-averaged model, in combination with a logarithmic velocity profile.Olsen and Skoglund [7] used a 3D k-ε turbulence model to predict the steady flow and sediment transport fluxes in a sand trap.A layer integrated 3D numerical model was developed to predict suspended sediment fluxes by Lin and Falconer [4].Olsen [8] and Rüther and Olsen [9][10][11] developed a full 3D model of flow and sediment transport with an unstructured grid to represent the domain.Their work focused on the formation of alternate bars and the initiation of meanders, starting from a straight channel.Different types of commercial numerical models, such as MIKE_3, TELEMAC-3D and Delft3D, have been used to predict the impact of hydraulic structures on the hydrodynamic, sediment transport and morphological processes.
Ruether et al. [12] used a 3D numerical model, named SSIIM (Sediment Simulation in Intakes with Multiblock Option), to predict the flow field and the sediment transport fluxes for the Kapunga water intake, in Tanzania.The model was used to predict the suspended sediment distributions in the flow approaching a water intake.Khosronejad et al. [13] developed a full 3D hydrodynamic model, using a finite volume method to solve the Reynolds-averaged Navier-Stokes (RANS) equations, and combined it with a 3D sediment transport model.Souza et al. [14] compared laboratory sedimentation experiments for a shallow reservoir against predictions obtained using a 2D numerical model, with the depth-averaged Navier-Stokes equations and a sediment transport algorithm.
For the numerical model simulations reported herein, predictions are made using a layer integrated 3D hydrodynamic and sediment transport model.In this model, a two-equation turbulence model is used, with details of the governing equations and boundary conditions being given.The original numerical model was developed by Hakimzadeh and Falconer [15] for simulating re-circulating flows in tidal basins.Refinements have been made to the model in order to predict more accurately the hydrodynamic and sediment transport processes in reservoirs.In the current study, the application of the enhanced numerical model is tested for simulating the scaled physical model of Hamidieh Reservoir.Numerical model results are compared to the data measured from the physical model, with the model predictions generally agreeing well with the measured data.

Governing Equations
The numerical model comprises two separate modules, including: a hydrodynamic module and a sediment transport module.The hydrodynamic governing equations are solved using a combined layer integrated and depth-integrated scheme [16], with the two-equation k-ε model being used for turbulence closure [15].The operator splitting technique and a highly accurate finite difference scheme are used to solve the suspended sediment transport equation [4].

Layer Integrated Hydrodynamic Model
The layer integrated numerical model used in this study was based on a model originally developed by Lin and Falconer [16].The governing mass and momentum equations were solved using the finite difference scheme on a regular square mesh in the horizontal plane and a layer integrated finite difference scheme on an irregular mesh in the vertical plane.
A sketch of the 3D grid and the relative positions of the governing variables in the x-z plane are illustrated in Figure 1.As illustrated, three layer types exist, including a top layer, bottom layer and middle layer.The top and bottom layer thicknesses vary with the x, y coordinates, to prescribe both the free surface and bottom topography, respectively.In contrast, the middle layers have a uniform (or non-uniform) thickness.The Equation of continuity in differential form for layer k is given as: Horizontal velocity

Vertical velocity
Turbulence model where wk−1/2 is the vertical velocity component at the interface between layer k and layer k−1, and the continuity equation is refined to give: where ζ denotes the water surface elevation relative to the datum.Similarity, vertical integration of the momentum equations for the k-th layer gives: where ± refers to the vertical (or z) elevations of the layer interfaces between the k+1, k and k−1 layers, respectively, = t time, are the layer-averaged velocities in the x, y and z directions, respectively and = ∆ , = ̅ ∆ are the discharges per unit width in the x and y directions, respectively.At the free surface (where = 1), the terms ( ) / and ( v ) / can be eliminated using the kinematic free surface condition and the Leibnitz rule.At the bed, the terms ( ) / and ( ̅ ) / become zero according to the no-slip boundary condition.

Turbulence Model
River flows are generally always turbulent, and the turbulent fluctuations contribute significantly to the transport of mass, momentum and heat.Even for the case of flow in a straight, smooth channel, the turbulence generated near the bed is highly complex, featuring a variety of three-dimensional coherent structures [17].The turbulence is complicated even further by geometrical variations, such as those due to bed forms, roughness elements and vegetation, changes in the river cross-section, bends causing secondary currents, confluences associated with strong shear layers and all types of man-made structures, such as dikes, bridge piers, groins, etc. [18].Strong turbulent flows usually occur in regulated reservoirs, primarily due to the high flow rates.
In the numerical model reported herein, the eddy viscosity concept has been used to represent the Reynolds stresses.The horizontal eddy viscosity was assumed to be constant along the water depth, with the depth-averaged eddy viscosity value being set for all layers.The depth-averaged horizontal eddy viscosity was calculated using [19], to give: where σ = the eddy viscosity coefficient, C = the Chezy roughness coefficient, H = the depth of flow and U, V = the depth average velocity components in the x, y directions, respectively.The vertical eddy viscosity was evaluated using the layer integrated form of the k-ε equations based on Hakimzadeh [20].The two equations were discretised using a semi-implicit scheme, and the resulting equations were then solved using the Thomas algorithm, as outlined in Hakimzadeh [20].
The following Equations were used for the free surface boundary conditions: and for ε , the following condition suggested by Krishnappan and Lau [21] was used: where the distance from the free surface to the grid centre, f k = the corresponding value for the turbulent energy, ε f c = 0.164, c μ = 0.09 and κ = 0.41, the von Karman constant.

Sediment Transport Model
Depending on the size and density of the bed material and the flow conditions, sediment particles can be transported in the form of bed load or suspended load.The 3D transport equation for the concentration of sediment in suspension is written as: where S = the suspended sediment concentration, Ws = the particle settling velocity and ε ,ε x y and z ε = the sediment mixing coefficients in the x, y, z directions, respectively.For suspended sand particles in the range of 100-1000 μm, the following equation was used to determine Ws [22]: where ν = the kinematic viscosity for clear water, S = the specific density of suspended sediment and Ds = the representative particle diameter of the suspended sediment particles.
The mixing coefficients were related to the turbulent eddy viscosity through the following Equations: where σ , σ h v = the Schmidt number in the horizontal and vertical directions, respectively, with values ranging from 0.5 to 1.0.
In solving Equation ( 8), four types of boundaries (inflow, outflow, free surface and sediment bed boundary) needed to be specified, with more details being given in Lin and Falconer [16].

Numerical Simulation
Details of the numerical schemes used in the hydrodynamic, turbulence and sediment transport modules can be found in Falconer et al. [23], Lin and Falconer [4,16] and Hakimzadeh and Falconer [15].
The transport of suspended sediment has been solved in the model by using a splitting algorithm, which splits the original 3D advective-diffusion equation into a 1D vertical equation and a 2D horizontal equation.An implicit finite volume method, used in the vertical sub-equation, was able to maintain stability for small layer thicknesses, and the ULTIMATE QUICKEST scheme was used to give highly accurate approximations of the horizontal advection process.For more details, see Lin and Falconer [4].

Project Background
The reservoir studied in this project is located 11 km away from Hamidieh Town, along the Karkhe River, downstream of the Karkhe Reservoir dam, in Iran (Figure 2).The value of the longitudinal slope of the Karkhe River in this area is approximately 0.00259.Currently, there are two water intakes, named the Vosaileh and Ghods intakes.The Vosaileh water intake channel is 10.8 km long, with a maximum discharge of 60 m 3 /s, while the Ghods intake operates through a 2.5 km-long channel, with a maximum discharge capacity of 13 m 3 /s (Figure 3).
Due to developments in the irrigation and drainage networks in the Azadegan and Chamran regions, the volume of water transferred from the existing intakes is no longer able to meet the water demand.Hence, it is necessary to increase the flow rate by building new water intake structures.
The planned Hamidieh dam is 192 m long, 4.5 m high and with 19 spillway bays and 10 sluice gates.A new water intake for Azadegan, with an inlet width of 56 m, eight opening bays and four under sluice gates, will be built to replace the Ghods water intake.The aim of the new scheme is to increase the discharge capacity from 13 m 3 /s to 75 m 3 /s.The Vosaileh water intake will also be replaced by a new water intake, with an inlet width of 86.6 m, 16 opening bays and 13 trash racks, to increase the discharge capacity to Chamran from 60 m 3 /s to 90 m 3 /s.Figure 3 shows a plan view of the Hamidieh regulated reservoir and the associated hydraulic structures.
An undistorted 1/20 scale physical model of the reservoir and its relevant structures was constructed, and experiments were undertaken to investigate the flow, sediment transport and morphodynamic processes in the reservoir.Hydrodynamic and sediment concentration data were collected to verify the model.The aim of the study was to better understand the operation of the water supply system and to improve the velocity structure and sediment transport processes in the vicinity of the structures.

Hydrodynamics
The numerical model simulations were conducted based on the laboratory measured data.Physical model experiments were undertaken for various operation scenarios of the intakes and sluice gates.The data included acquiring detailed velocity values within the reservoir and, in particular, in the vicinity of the hydraulic structures.Figure 4 shows the bathymetry of the reservoir and the intakes of the physical model.In the physical model study, the velocity was measured at several points in front of the intakes.The locations of these points are shown in Figures 5 and 6.The velocity measurements at these points, together with the observed water surface levels, were used for calibrating the hydrodynamic model.The velocity was measured using a current meter at 0.6D (D = water depth) below the water surface.This depth was chosen based on a logarithmic velocity profile assumption.In applying the 3D layer integrated model, the water depth was divided into five layers.The thickness of each layer was set to 0.057 m, except that the top layer thickness was initially set to 0.051 m.The horizontal grid size was 0.25 m.The water surface level was calibrated by comparing the numerical model-predicted water surface levels to the observed levels.Table 2 lists the velocity values obtained from both the physical and the numerical models for the locations shown in Figures 5 and 6.Table 2 shows that the model-predicted velocities generally agree well with the measurements for this complicated design.The average error was about 11%.The results obtained from both the physical and numerical models are presented for the undistorted 1:20 scale model of Hamidieh Reservoir.The flow patterns and velocity contours in the reservoir, predicted using the numerical model, are illustrated in Figure 7     Var.
From Figure 8, it can be seen that the numerical model-predicted velocity values agree reasonably well with the measured values, although some variations can be observed.The type of turbulence model has an important role in the accuracy of the predictions.
After calibrating the numerical model, using measurements from the physical model (Scenario S1), two more scenarios were considered to acquire a better understanding of the flow and velocity distributions in the regulated reservoir.The key characteristics of these scenarios are presented in Table 3.The numerical model is switched to one layer for these two scenarios (i.e., depth-integrated).The numerical model-predicted flow pattern and velocity contours in the regulated reservoir for these two scenarios are illustrated in Figure 9a,b.As can be seen from the results (Scenarios S2 and S3, Figure 9a,b), when the sluice gates can release more water than the flow velocity near Chamran, intake rises significantly and the performance of this intake is improved.Simultaneously, the possibility of sedimentation in the vicinity of this intake is reduced.The numerical model results reveal that variations in the flow velocity near Chamran intake, due to an increasing reservoir discharge, are more sensitive at Chamran in comparison to such changes for the Azadegan intake; with this thought to be due to the lower water depth and bathymetry in front of Chamran intake.
An analysis of the hydrodynamic results identifies three main issues: -A non-uniform velocity distribution zone (in the horizontal plane) exists on the left-hand side of the bank upstream of Chamran intake.This is also observed in the physical model (see Table 2, A-M and E-N points; and see Zone A in Figure 7).-The velocity pattern for the regulatory condition over the artificial dredging zone (in front of Azadegan intake) indicates that this is a location of high potential sediment deposition (see Zone A in Figure 9a and Zone B in Figure 9b).-The flow velocity field near the Chamran intake shows more variations in comparison to the velocity field near Azadegan intake (see Figures 8 and 9a,b).

Sediment Transport
In order to model the sediment transport processes, the type of sediments (i.e., cohesive or non-cohesive), the sediment particle size distribution and an understanding of the dominant processes of sediment transport (i.e., suspended load or bed load) are required.Previous studies have shown that suspended sediment transport is the dominant process for this case [24], and the majority of the sediments found in this region belong to the non-cohesive type.
Details are given below of the three sediment transport scenarios considered in this study: 1. Low flow without sluice gates (S4).
3. High flow with sluice gates along both sides of the dam in operation (S6).
The attributes of the scenarios are presented in Table 4.
Measurements of sediment concentrations were made at five sites, with the locations of these sites being shown in Figure 10.The samples of suspended sediment concentrations (SSCs) in the physical model were taken at several points along the water depth and then depth averaged.Numerical model simulations were undertaken to predict the suspended sediment concentration distributions for these three scenarios.Table 5 and Figure 11 show comparisons between the measured and calculated SSCs at the locations shown in Figure 10.The model-predicted patterns of the suspended sediment concentrations are shown in Figures 12-14, for Scenarios S4, S5 and S6, respectively.The SSCs computed using the numerical model agreed well with the concentrations measured in the physical model.The maximum error was less than 21.3%, and the average error was 8.28%, both of which were considered acceptable when considering the overall complexity of the sediment transport processes in the reservoir.The movement of suspended sediments in the reservoir, over the intakes and sluice gates, can be clearly seen in Figures 12-14.The results show that the suspended load has an important role in the sediment transport processes and is closely related to the particle size.Good agreement has been observed between the SSCs obtained from the numerical and physical models.The analysis of the sediment transport results identifies that: -The SSC in the Chamran intake is higher than that in the Azadegan intake (see CCH and CAZ in Table 5).This is due to the fact that the Chamran intake was located along the convex bank, where the flow velocity was greater than that at other locations across the cross-section.-The pattern of the suspended sediment movement in the reservoir shows that the sediment flux rate entering into the Chamran intake is not significantly affected by the intake gate and sluice gate operation schemes.This means that the sluice gates on the left-hand side do not have a considerable effect on reducing the SSCs passing through Chamran intake.-The problem of sediments entering the intakes (especially the Chamran intake) still exists.The performance of these hydraulic structures is only part of the reason for this problem.The dominant natural process of sediment transport (especially for the suspended load part) and the size of the sediment particles play an important role in connection with this problem.In mitigating against this problem, then, dredging and sediment control measures should be periodically undertaken at the upstream reach of the dam.

Conclusions
Details are given of the refinement and application of a three-dimensional numerical model to predict the flow and sediment transport processes through a reservoir.The aim of the study was to obtain an improved understanding of the ability of such a numerical model to simulate these complex processes for a real case study.An enhanced turbulence model was used to represent the complex flow patterns in the vicinity of the hydraulic structures.The corresponding numerical model predictions were compared to the scaled physical laboratory model data.The calibrated numerical model used in this case study could provide a platform for investigating the hydrodynamic and sediment transport processes in similar practical engineering studies.
A series of scenario runs of the operation schemes for the regulated Hamidieh Reservoir have been undertaken to investigate the governing hydro-morphological phenomena occurring in the system.The results obtained from the numerical model simulations suggest that: i.
The 3D layer integrated scheme used in the numerical model is capable of representing the governing three-dimensional hydrodynamic and solute transport processes with an encouraging level of accuracy for the cases considered.Due to the important role of the near-bed layer (with high sediment concentrations), the sediment transport processes obtained using the 3D layer integrated model is applicable for predicting the hydro-morphological processes in the natural environment.By using this method, the computational time has been reduced significantly.ii.
The computed SSC distributions predicted using the numerical model for Hamidieh Reservoir agreed well with the measurements taken from the physical model, with an average error of less than 8.3% (see Table 5).To acquire this level of accuracy, the hydrodynamic computational scheme and the selected turbulence model were found to be crucial.The results of this research showed that the selection of the k-ε turbulence model in reservoirs, where the Reynolds number is high, but without any significant local vortex near any hydraulic structures (such as gates), was the correct choice.An accurate estimation of the eddy viscosity had a direct impact on the calculated flow velocity distribution (both horizontally and vertically), and the mixing coefficient was found to be a key parameter in calculating the sediment transport fluxes.iii.
The calibrated numerical model, set up for Hamidieh Reservoir, revealed many hydraulic aspects that were generally in good agreement with the physical model study results.Some of the important aspects are summarised below: -A non-uniform velocity distribution zone existed upstream of Chamran intake on the left-hand bank side.The flow velocity field near this intake showed more variations in comparison to the velocity distribution for Azadegan intake.-The recommended artificial dredging zone, located in front of the Azadegan intake, did not affect the hydraulic behaviour significantly.-The location of Chamran intake had a considerable impact on the SSCs moving through the intake.-The intake and sluice gate operation schemes had a relatively small impact on the suspended sediment fluxes entering Chamran intake.
iv.The results of this research have indicated that the calibrated numerical model could be used in the hydro-morphological design, operation and management of regulated reservoirs, particularly following the application of the model to the study reported herein.By using this approach, laboratory tests were reduced, and thus, time and budget considerations were optimized.

Figure 1 .
Figure 1.Coordinate system for layer integrated equations.

Table 1
lists the hydraulic parameters used to set up the numerical model.These parameters were based on data obtained from the physical model experiments.

Figure 4 .
Figure 4. Bathymetry of Hamidieh Reservoir, at the physical model scale.

Figure 5 .
Figure 5. Location of velocity measuring points in the vicinity of the Azadegan intake.

Figure 6 .
Figure 6.Location of velocity measuring points in the vicinity of the Chamran intake.
(for Scenario S1).The speed profiles along the depth predicted by the numerical model for Scenario S1, at four points near the intakes (B and F for Azadegan and M and N for Chamran), are shown in Figure 8.The measured velocity values at 0.6 below the water surface in the physical model are also shown in this figure.

Figure 7 .
Figure 7. Numerical model-predicted flow pattern and speed contours for Scenario S1.

Figure 8 .
Figure 8. Model-predicted velocity distributions and measured velocities at 0.6D below the water surface for Scenario S1 (Var = variation, Az = Azadegan intake, Ch = Chamran intake, Nu = numerical, Ph = physical and TD = total depth (i.e., calculated by 3D layer integrated numerical model along the whole depth)).

Figure 9 .
Figure 9. Numerical model-predicted velocity pattern and speed contours for (a) Scenario S2 and (b) Scenario S3.

Figure 10 .
Figure 10.Location of sediment concentration sampling points (CAZ = concentration in Azadegan intake, CCH = concentration in Chamran intake and CRR, CRM and CRL = concentrations at right, middle and left sampling points).

Table 2 .
Comparison between model-predicted and measured velocity values for Scenario S1.

Table 3 .
The characteristics for extra hydrodynamic scenarios.

Table 5 .
Suspended sediment concentrations (SSCs) measured and predicted at the sampling points.