Mathematical modelling and interpretation of results for siltation phenomena of rivers in the vicinity of reservoirs

. Worldwide the processes of erosion, transport and sedimentation have a significant social, economic and environmental impact. In the past decade the most important external factor that has significantly increased the erosion and sediment process is the human factor throughout its daily activities. The changes of land cover and use, the intensive deforestation, urbanization and infrastructure development, sand extraction from river beds alter the natural watercourse and thus leading to erosion, sediment transportation and finally to the silting of rivers and reservoirs. The process of erosion and sedimentation occurs, among other areas, in the hydro-technical constructions with significant effects. In the next 30 years it’s expected to lose more than 50% of the storage capacity of the reservoirs in the world due to the silting phenomenon. The control and monitor of the sediment dynamics represents a critical component in the sustainable management of water resources. Also the understanding of the involved phenomena is essential in order to develop solution for reducing the sediment and erosion that occur and adopt measures to extend the life of the reservoirs. This paper taps the solution of mathematical modelling using dedicated software in order to establish the morphological changes that occur in Ogrezeni-Crivina reservoir, located on the river Arges.


Introduction
Sediments represent solid materials in the form of particles which after being detached from the earth's crust (river bed) through the action of natural forces and the transport in streams are deposited at the surface of the earth. The sediment process consists in some major steps: the erosion process which consists in the dislocation of the soil, in particles under the action of water and other agents; the transport process with the water acting as a movement force; the sedimentation process of the transported particles.
The erosion process is generated by the dislocation and transport of solid elements in the soil, mainly particles, through waters, glaciers, wind and waves. A particle subjected to this process tends to move as the result of the action of an active force and it opposes by its weight and friction forces. In the case of water the hydrodynamic forces represent the active forces and the opposing forces are defined by the frictional force corresponding to the underwater weight. Generally in mountainous and hilly areas, erosion is predominant, while sedimentation occurs almost all of the time in the plain areas.
In order to analyze the sediment transport process in the water bodies, it is considered that the material is noncohesive whose particles move practically independent, with no phenomena occurring between them that can lead to the adhesion of the particles. Another hypothesis that is considered is that the flow in rivers is turbulent. The sediment movement into a stream of water takes place in jumps. In a turbulent flow upon a particle of sediments from the riverbed drive both longitudinal and ascension forces. The longitudinal forces depend on the particle strength (shape / friction) and the vertical velocity distribution in the water course. The ascension forces are the result of the pressure forces that act on the lower and upper faces of the particles. Resistance forces that oppose the movement come from the particle weight under water and from the interaction with the nearby particles. In the terms of particle movement we can distinguish the following: sliding movement, rolling movement, leaping movement and floating movement. Also the major characteristics of a particle are its dimension, from, specific weight and velocity of drop.

Case study
The case study focuses on the reservoir on Arges River created by dam Ogrezeni-Crivina with the main purpose of supplying drinking water to Bucharest city, the residence of approximately 2 million inhabitants. The complex construction from Ogrezeni-Crivina is composed of a main dam equipped with movable weirs and flushing gates and also a water intake with the purpose of taking water from Arges River and transporting it gravitationally to Bucharest.
Since 1905 Arges River is one of the primary Bucharest water supply sources. Consequently, it was initially planned to supply to Bucharest more than 700.00 m3 per day, from which approximately 150.00 m3 per day would be treated for drinking water and the remaing for industrial use, using the surplus to improve the water quality of Colentina River and dilution of Dâmboviţa River, both through diversion works.
The works of the complex layout began in the year 1939 with the open channel diversion, finished in 1943 for approximately 13.50 km when construction was interrupted. Channel works resumed in the year 1946 and it was put into operation in 1950, being used only for industrial water supply, the main consumer being the Grozavesti Electric Power Plant -thermoelectric.
At the same time the construction of the dam on the Argeş River between the villages of Crivina and Ogrezeni began. On the left bank intake, four sand basins were implemented, with a total capacity of 10.00 m 3 /s in order to supply the capital. The development of the city and implicitly the increase of water demand led to the expansion of the water intake and to the construction of another 6 longitudinal sand basins with a capacity of 15.00 m 3 /s and two pumping stations for supplying a flow of 3.00 m 3 /s to the water treatment plant in Arcuda.
Today Argeș River, in the considered section of the Crivina-Ogrezeni dam, can provide a total flow of approximately 18.80 m 3 /s distributed as followed: • 40% for drinking water demand; • 70% for industrial water demand; • 95% dilution flow for Dambovita River; • 70% for recharging the cascade lakes on Colentina River. The movable weir was designed with 6 openings ( Fig. 1) of 26.00 m each, equipped with segment gates, which are operated by means of drive mechanisms, located in the cabinets inside the concrete piles of the dam. It also includes 3 washing openings of 20.00 m length with 6 flat surfaces in order to flush the channel in front of the intake. Each gate is equipped at the top with a lower flap to evacuate the floating bodies, snow and ice during the winter. The dam has the following characteristic [1]: A. Characteristic volumes: 104.00 mASL. In order to dissipate the hydraulic energy of water the stilling basin, with a length of 13.00 m, was built downstream of the dam and consisted in a concrete slab and wooden beams. Because the wooden beams were destroyed, a second dissipation compartment was built. The basin was designed with reinforced concrete teeth with square section, placed in the shape of chess board.
Over time, however, studies have shown the need to rehabilitate the dam and reservoir. A project was drawn up and, starting with 1981, the execution began and the works ended in 1984. The new dissipater includes: • Reinforcements with cement injections; • The second compartment provided with prismatic teeth located on a single row; • The hydraulic jump limitation beam; • River boulders downstream the stilling basin.
The exploitation of the dam over the next 10 to 15 years has highlighted the fact that the downward scouring continued on both the general and the adjacent area of the construction, reaching about 4,00 -5,00 m.
Measurements and recordings regarding the total suspended solids have been carried out in the Arges River since 1940, mainly due to its vital role for water supply. To determine the values of the sediment present in the river (dragged, suspended and total amount), the historical recordings performed at the hydrometric stations were analyzed in order to establish a correlation, for the considered area, between the solid flow and the liquid flow.
The backwater calculations based on the topobathymetric profiles revealed the following: • The degree of sedimentation in the reservoir is significant, therefore not allowing flood attenuation (only the transit of the flood). Currently the lake is silted in proportion of 80%, with a water depth of approximately 2 meters just near the dam and intake; • During a flood event the retention structure can barely provide the discharge of the verification flow; • The amount of sediments deposited on the right bank has the appearance of a peninsula with a level close to the Normal Pool Level and redirects the flow to the intake, thus favoring the operation conditions, especially during winter periods. • During flood conditions, due to the massive sedimentation of the reservoir, the water levels increases significantly, 125 m upstream the dam, reaching approx. 0.5 m for flows up to 1000 m3/s and approx. 1 m for verification discharge. In this case, for flows over 1000 m3/s, Crivina dam is not operated in safety conditions, the water level exceeding the level of the banks from both sides, with the possibility of flooding the adjacent areas. The sedimentation in the upstream area of the reservoir has been favored by the operation regime of the dam, without adequate maneuvers of the weirs that could have flushed the sediments and maintain the storage volumes to preserve the hydraulic equipment in good working conditions. In addition, the preparing maneuvers for winter periods amplify the sedimentation process due to the reduction of the rivers upstream section in order to increase the velocities to facilitate the water transit and to avoid the blockage of the intake.
Today, the deposits of solid materials prevent the free transit of the influent flow into the central area and the right bank of the dam.

Mathematical modelling
There various mathematical methods developed to simulate the sedimentation process inside reservoirs based on the present advance computational technology and processing capacity. All the software models include three basic components: flow, sediment transportation and dedicated functions. The mathematical models used to study the reservoir erosion, transport and sedimentation processes are cheaper, less time consuming, user friendly interface and easier to manage than the physical modelling and especially than the empirical methods. In our particular case, the mathematical model is conceived using DELFT3D in order to analyze the river's morphology.
The DELFT3D software consists of various integrated modules that together allow hydrodynamic flow simulation (assumption of shallow water), the salinity and the heat of the water, sediment transport and morphological changes (erosion, silting) as well as the modeling of ecological processes and water quality parameters. It was developed in a multi-dimensional hydrodynamic simulation program (2D and 3D), that determines the flow and transport processes that result from tidal and weather action.
The functionality of the program takes into consideration the momentum equation, the continuity equation and the transport equation.

Continuity equations
The depth-averaged continuity equation is given by [2]: S represents the contributions per unit of area due to the water withdrawal or discharge as well from evaporation or perspiration.

Momentum equation
The depth-averaged equations in the case of momentum preservation in X and Y directions are given by the following relations [2]: On X direction: On Y direction: Terms Px and Py take the following form: Mx and My come from external sources such as hydraulic forces, water discharge or withdrawal, wave effects.

Sediment transport equation
The sediment transport is determined by solving the three-dimensional advection-diffusion equation [2]: Local velocity magnitudes are based on the results of hydrodynamic calculations and for sediment fractions, a mixture of "sand" and "mud" can be selected. The erosion and sedimentation rates of each particle type are calculated. The simultaneous presence of several sediment fractions has implications both for calculating the local sedimentation rate and for calculating the density of the mixture.
Each particle has its own settling rate depending on space and time. Transport of suspended solids is calculated for all sediment fractions (Van Rijn 1993). For waveless simulations, the magnitude of the total load on a horizontal rivel bed is deterimined using a formula given by Van Rijn (2000) [2]: The net change of sediments due to transport in suspension is calculated as follows:

Sediment exchange with the river bed
Van Rijn's (1993) approach is used for sediment particle fractions, and reference heights are determined based on the bed roughness. The sediment concentration at the reference height is determined using an adapted formula of Van Rijn (1984) to include the presence of multiple sediment fractions [2].

Model design
The bases of the mathematical model represent the topobathymetric and hydrometric measurements that were made in order to achieve the study. As any modelling software DELFT3D needs 3 basic things: geometry (meshing, bathymetry), boundary conditions (inflow, outflow, water depth / level etc.) and observation / control points in order to calibrate the model and compare the obtained results. The grid was constructed as a 2D domain mesh and afterwards the level of the river bed (from bathymetry measurements) was assigned to the each node to simulate the actual state of the sediments. The messing process was done with curvelinear elements and was optimized using triangular interpolation, smoothing functions in order to achieve the best ration between run time and result precision. In this particular case a small cell size can lead to a run time of approximately 45 -60 days and a large cell size can cause anomalies and errors in the model. The boundary conditions consist in the inflow from Arges River (hydrometric records) and at the outlet the imposed level from the operating conditions of the movable weir and also the intake from the water treatment plant. The sediment yield was provided by the reservoir operator and also by studied conducted in the same area [1,3].

Fig. 4. Morphological grid with assigned water depths and boundary conditions (mASL).
Due to the fact that the only interest was to assess the evolution of the sediments (erosion and deposits) and the preferential flow paths, the limits of the rivers (right and left bank) were modeled as wall conditions so that the total flow will be transited thought the river bed without adjacent flooding.
The simulation run time of the model was 365 day, more precise for the entire year 2016, for which the needed data was of very good quality in the meaning hydrometric measurements and specific sediment load. The alternation of dry weather flows with very low inflows produced the sedimentation process and also rain events with peak flows produce the erosion and sediment transportation process (Fig. 5). Also the specific load of sediment was very high in this specific section of the river, a total amount of suspended solids of 25,6 kg/s, value that significantly increase during rain events ( Table 1). As mentioned above an improper cell size can cause instabilities in the model. An important fact that was taken into consideration was the simulation time step which can cause instabilities of the flow. In consequence, based on Reynolds and Courant numbers the time step of the simulation was calculated at a value of 1 second and furthermore to improve the accuracy of the results the morphological update time was considered less than a minute. Given the size of the model 1600m x 70m, the small meshing and the 1 second time can lead to time consuming simulation and one trial can last almost 7 days.

Results
As in can be observed in Fig. 6 -2 months simulation period and Fig. 7 a year simulation period the changes of the river bed morphology are significant. After a couple of months it can be observed that the preferential flow tends to amplify the meander belt of the river due to the erosion process that occurs. Also the sediments are deposited in the concave area of the river due the orientation of the velocity vectors. The total amount of sediment mass that resulted after the simulation was approximately 126 000 tons/year. In these circumstances the total volume of sediment deposited over a year in the reservoir is 60 000 mc/year, volume which is lost from the available storage capacity. The velocity maps reveal that during dry weather the reduced flow causes low velocities in the river, at an average value of 0.3 m/s, favoring the sedimentation process (Fig. 8). On the other hand, rain event induce, during peak flows, a significant increase of the velocity magnitude thus favoring the transport and erosion process. Furthermore in narrow cross sections the velocity reaches values up to 2.6 m/s causing the erosion of the left and right bank and also the river bed (Fig 9). In consequence it can be observed that the velocity vectors and magnitude shape the river in it's entirety from the river bed to the river banks.   Figure 9. Average depth-averaged velocity for peak flows (m/s).

Conclusions
From the simulation it can be concluded that the reservoir losses almost 60 000 mc/year from the storage capacity. During dry weather a significant amount of deposits are formed in the river, building a preferential flow path, and during rain events with peak flows the morphology suffers important changes due to the erosion and transport processes that occur. Furthermore the analyze reveals that Arges River is subjected to heavy sedimentation and without proper maintenance and operation it can lead to the silting of the Ogrezeni-Crivina reservoir, endangering the drinking water supply of Bucharest. The obtained results can be used to for the development of solutions such as diversion works, flushing maneuvers, removal plan for the sediment deposits etc. During flood events, with peak flows, when the velocity increases significantly the weirs can be