Three-dimensional numerical modeling of hydraulics and morphodynamics of the Schwarzenbach reservoir

In this study the numerical model SSIIM 2 is used to simulate the flow situation and the suspended sediment transport within the Schwarzenbach reservoir in Germany. Hydrodynamic simulations are carried out to assess the influence of wind forces and different discretization schemes on the calculated flow field. A hydraulic plausibility check is performed based on stationary ADCP measurements to assess the performance of the model. Both the consideration of the wind and the plausibility check using an ADCP are hardly used in large reservoirs so far. The simulation results show a complex flow field with two large (re)circulation zones in the middle of the reservoir, whereby the temporal development of the simulated and measured velocities have comparable characteristics. Moreover, morphodynamic simulations are performed to compute the suspended sediment transport. The results show that the settling behavior of the sediments is mainly influenced by the inflow discharge and the operation level of the reservoir.


Introduction
Three-dimensional numerical modeling has become increasingly important in engineering during the last decades.One reason is that the application of one-or two-dimensional numerical models is for many practical cases not satisfying due to complex occurring flow situations.Especially in morphodynamic modeling it is essential to take these threedimensional effects into account to ensure correct model results.In this study the numerical model SSIIM 2 is used to simulate the flow situation, the suspended sediment transport and the morphological processes within the Schwarzenbach reservoir in Germany.Several algorithms and parameters, such as different upwind schemes and the influence of wind, are tested to ensure correct hydrodynamic boundaries for the morphological simulations.The simulated results are validated with stationary ADCP measurements, conducted during different reservoir operation conditions throughout the year 2016.Finally, different hydrodynamic scenarios, with respect to inflow and operation level, are used to gain knowledge regarding the sediment transport within the reservoir for different operation conditions.

Site description
The Schwarzenbach reservoir, located in the northern Black Forest, Germany, was built in 1926 to serve as upper reservoir in Europe's first large scale hydro-electric pump storage system.The reservoir has a length of 2.2 km, a maximum width of 600 m and a maximum depth of 47 m, at maximum operation level, leading to a total storage capacity of 14.4 million m³.The permissible operation level is between 628 m.a.s.l and 668.5 m.a.s.l [1].
Figure 1 shows the measured bed elevations (2012), the shore line at maximum operation level, the operation outlet and the three tributaries, which feed the reservoir from two catchments with a total area of 50 km² [1].The inflow to the Schwarzenbach reservoir consists of a pumped water volume and a natural inflow volume.The two creeks, Schwarzenbach and Seebach, enter into the impounded reservoir at its head.The Raumuenzachueberleitung (RUB) is a transition tunnel that collects several headwater streams from the adjacent catchment.As shown in Figure 1 the Raumuenzachueberleitung mouths laterally into the reservoir.The average annual inflow of the three tributaries is approximately 2.1 m³/s.With 55 % of the total natural inflow discharge, the Raumuenzachueberleitung has the largest share of the inflowing water volume, whereas the share of the Schwarzenbach (35 %) and the Seebach (10 %) is rather small [2].Due to the pluvio-nival discharge regime of the tributaries, the inflowing discharge during summer is significantly lower than in winter and spring.Additionally to the natural inflows, water from the river Murg is pumped into the reservoir.The percentage of the pumped water of the total inflowing water quantity is 36 % (2005-2015; [2]), inducing strong water level fluctuations due to pump operation.However, the average natural inflow is higher than the amount of pumped water, which underlines the hydraulic peculiarity.

Numerical method
For the numerical investigations the fully three-dimensional numerical model SSIIM 2 is used [3].SSIIM is an abbreviation for 'Sediment Simulation In Intakes with Multiblock option' and has been successfully applied in previous studies to model the flow field, suspended sediment transport and depositions in reservoirs [4][5][6][7].In order to compute the water flow, the CFD-program solves the Reynolds-averaged Navier-Stokes equations (RANS) in three dimensions.The finite volume method is used as discretization scheme to transform the partial differential equations into algebraic equations.An implicit time discretization is implemented in SSIIM so that the use of large time steps is possible.The wind induced shear stresses (  ), acting on the surface of the water body, are calculated as follows: where ρ a is the mass density of air and U a,10 is the sustained wind speed 10 m above the water surface.The friction coefficient (  ) in this study is calculated according to van Dorn [8].Thus, the friction coefficient is 0.01 for wind speeds less than 5.6 m/s and increases at higher wind speeds.
The model uses an unstructured, non-orthogonal and adaptive grid, which is made out of tetrahedral and hexahedral cells.The implemented algorithm for wetting and drying allows the grid to move with changes in the bed and water levels.Consequently, only the water phase is modelled, and the number of cells might change after each time step.For maximum operation level, the generated grid consists of 81,300 cells in total (spatial resolution 10 m x 8 m).Up to 18 vertical cells are generated in the deepest areas of the reservoir.The variables at the three inflow boundaries are defined by Dirichlet boundary conditions, whereas a zero-gradient boundary condition is used for the intake.The velocity profiles near the bed are described by the wall law for rough walls [9].
The computation of sediment transport is divided into suspended and bed load transport.The suspended sediment transport is calculated by solving the transient convectiondiffusion equation.If different sediment sizes are implemented, the equation will be solved for each sediment size separately [3].Furthermore, the bedload transport is calculated by an empirical formula according to van Rijn [10].The roughness is chosen with 0.02 m and the active layer thickness to be d max .

Hydrodynamic simulations
Hydrodynamic simulations are carried out in a first step to assess the influence of different algorithms and parameters on the simulated flow field.Among those are simulations to investigate the effect of wind on the flow field as well as the effect of the Power Law scheme (POW) compared to the Second Order Upwind scheme (SOU), which are used for discretization of the convective terms of the Navier-Stokes equations.The simulations are performed for steady state conditions.The chosen conditions are: the power plant is in turbine operation (16.8 m³/s) and the natural inflows are slightly below the annual maximum (28.2 m³/s).The water level is kept constant on a level of 680.0 m.a.s.l.Due to the small simulation period and the large initial volume discontinuity is neglected.
To assess the performance of the model, a plausibility check is performed based on stationary ADCP measurements from the year 2016, which were performed by the Limnological Institute of the University of Konstanz within the CHARM project [11].The ADCP was installed in an upward-looking configuration on the bottom of the reservoir approximately 270 m northwest from the operation outlet (see Figure 1).Two unsteady periods including turbine and pumping operation are simulated with high temporal resolution.

Morphodynamic simulations
Morphodynamic simulations are carried out to compute the suspended sediment distribution within the reservoir.Different scenarios are simulated, which differ in terms of operation levels and inflow conditions.The power plant is in turbine operation, with the turbine flow rate being equal to the reservoir inflow.Five different grain fractions are implemented into the numerical model.The grain diameters correspond to the averaged percentile values d 10 , d 30 , d 50 , d 70 and d 90 of different sediment samples taken from the reservoir [12].To calculate the fall velocities of particles smaller than 0.1 mm Stokes law was applied [13].The fall velocities of larger particles are calculated according to Ahrens [14].Due to an upstream located sediment trap, the transition tunnel (RUB) does not carry any sediment into the reservoir.Sediment input occurs only through the tributaries of the Seebach and the Schwarzenbach.For the simulations the mass concentration for each grain fraction is set to be 10.6 mg/l.Since no information regarding the sediment concentration in the inflows is available, the concentrations and the simulation period are chosen in a way that the tributaries do not dry out during the simulation due to morphological bed changes.

Hydrodynamic simulations
All hydrodynamic simulation results show a complex flow field with two large recirculation zones in the middle of the reservoir, mainly as a result of the lateral inflow from the transition tunnel (RUB).As it can be seen from Figure 2, the highest depthaveraged flow velocities occur in close vicinity to the RUB and the Schwarzenbach.The simulated flow field using the SOU scheme differs in individual areas of the reservoir significantly from the solution calculated by means of the POW scheme.Figure 2 shows that the flow resulting from the Schwarzenbach is more concentrated and therefore higher flow velocities occur when applying the higher order scheme (SOU).The reason for the differences is mainly the false diffusion occurring when applying the POW scheme.The false diffusion is particularly large in the north-western area of the reservoir because the main flow is oblique to the grid.In addition, when using higher-order discretization schemes, vortices and recirculation zones are more pronounced or sometimes additional eddies occur.As shown in Figure 2, when applying a higher-order discretization scheme, two recirculation zones become visible in the north-western area of the reservoir.
Additional tests were performed with the QUICK and the Van Leer upwind schemes, but only marginal differences compared to the SOU scheme could be seen.
According to a measurement campaign and the nearby weather station Hornsigrinde (48.61°N, 8.20°E), the prevailing wind direction in the area of the Schwarzenbach reservoir is West-South-West (247.5°).The wind direction as well as the average wind speed (1.57m/s) is implemented in the model and remains constant over time.The investigations demonstrate that the wind influences the flow field within the reservoir under the tested hydraulic boundary conditions.It can be seen from Figure 3 that the flow field within the reservoir is divided into two eddies as result of the inflow of the transition tunnel.The flow, coming from the RUB is split at the opposite bank and moves partly in northwest and partly in southeast direction.The flow in northwest direction runs initially along the north shore.Without the influence of the wind, the flow detaches from the shore before reaching the bay.In addition, a circulation zone is formed across the entire width of the reservoir.The consideration of wind in contrast prevents the flow from detaching and thus prevents the formation of the circulation zone.By comparing surface and bed near velocity vectors, it can be seen that wind has an influence on the flow field down until the bottom, even if the largest changes of the flow field occur surface near.
Comparisons between the measured and simulated flow velocities illustrate comparable patterns and characteristics.Figure 4 shows the measured (ADCP) and the simulated flow velocities (depth averaged) in direction north for a period of 6 days.The presented results are simulated using the POW scheme without considering the wind.Amplitudes with a positive velocity represent pumping operation, whereas negative values represent turbine operation.The model's response to changing hydraulic boundary conditions is comprehensible and tends to be well aligned with the measured velocities.The northwards flow velocity increases during pumping operation (Figure 4) and the eastwards velocity increases during turbine operation in both, the simulation and the measurements.However, the numerical model underestimates the velocity magnitude during pumping and turbine operation.Moreover, a consideration of the wind and the application of a higherorder discretization scheme do not improve the correlation between the measured and simulated velocities.Possible reasons may be that the model neglects temperature-related density effects and that it overestimates the fluid viscosity.In addition, the measurement of low occurring flow velocities is prone to errors and could also account for the discrepancies.

Morphodynamic simulations
The results of the morphodynamic simulations demonstrate that the fine particles are carried further into the reservoir, while the coarse particles deposit close to the inflow area (Figure 5).Thus, there is an increasing refinement of the sediments into flow direction.Figure 5 shows the suspended sediment concentrations close to the surface for the case of a maximum annual inflow at full operation level, whereas Figure 5b shows only fraction 5 (d = 5.1 µm).This decreasing grain size gradient towards the dam is characteristic for reservoir sedimentation and corresponds to both, the theoretical expectations [15] and the bed sediment samples taken from the reservoir [12].6 (mean annual inflow and maximum annual inflow; mean as well as maximum operation level).The results illustrate that the suspended sediment distribution is mainly influenced by the total inflow discharge as well as by the operation water level.If the inflows correspond to the mean annual maxima, larger flow velocities occur in the entire reservoir and the particles are transported further into the reservoir.Moreover, there are fewer sediments in suspension at lower inflow discharges.Figure 6 shows also that the suspended sediment distribution is influenced by the water level within the reservoir, but with a slightly smaller effect compared to the inflow discharge.This is mainly due to the smaller flow cross-section and the resulting higher flow velocities at lower operating levels.Furthermore, it can be seen that the water level (and the inflow discharge) affects the shape of the sediment plume and thus the resulting deposition pattern.In contrast to the natural inflow discharges, the flow coming from the RUB and the water level is influenced by the operation of the dam.Since the system is sensitive to changes in water level, the operator can (to some extent) control/manage the deposition pattern of the sediments.For example, the sediment plume runs along the south shore at full operation level, whereas the concentration of suspended sediment at lower water level is distributed relatively evenly over the cross-section in the northwestern area.

Conclusion
Within this study the application of a 3D numerical model to simulate the flow situation and the suspended sediment transport in the Schwarzenbach reservoir is presented and discussed.It can be seen, that a complex flow field exists in the reservoir due to the lateral inflow of a transition tunnel.A closer look to numerical algorithms shows, that the wind influences the flow field to a high degree in the whole reservoir, whereas the application of a higher-order discretization scheme leads to different solutions only in individual areas.The model is verified by means of stationary ADCP measurements.This method has due to the small flow velocities in large reservoirs hardly been used so far.Moreover, the morphodynamic simulations demonstrate that the suspended sediment distribution is primary influenced by the inflow discharge and the operation level of the reservoir.Consequently, by controlling the water level it is possible to influence the deposition pattern within the reservoir.

Fig. 2 .
Fig. 2. Comparison of the simulated flow fields (depth averaged velocities) by applying the POW scheme (left) and the SOU scheme (right).

Fig. 3 .
Fig.3.Comparison of the flow field (POW) without (left) and with (right) consideration of the wind forces.

Fig. 4 .
Fig. 4. Comparison of the simulated and measured flow velocities towards the north.

Fig. 5 .
Fig. 5.Total suspended sediment distribution (left) and suspended sediment distribution of grain fraction 5 (right) in the top cell

Fig. 6 :
Fig. 6: Suspended sediment distribution of grain fraction 5 depending on the operation level and the inflow discharges

Table 1 .
Fall velocities and characteristic diameters of the individual grain fractions.