3D Numerical Simulation of Gravity-Driven Motion of Fine-Grained Sediment Deposits in Large Reservoirs

: Deposits in dam areas of large reservoirs, which are commonly composed of ﬁne-grained sediment, are important for reservoir operation. Since the impoundment of the Three Gorges Reservoir (TGR), the sedimentation pattern in the dam area has been unexpected. An integrated dynamic model for ﬁne-grained sediment, which consists of both sediment transport with water ﬂow and gravity-driven ﬂuid mud at the bottom, was proposed. The incipient motion driven by gravity in the form of ﬂuid mud was determined by the critical slope. Shallow ﬂow equations were simpliﬁed to simulate the gravity-driven mass transport. The gravity-driven ﬂow model was combined with a 3D Reynolds-averaged water ﬂow and sediment transport model. Solution routines were developed for both models, which were then used to simulate the integral movement of the ﬁne-grained sediment. The simulated sedimentation pattern agreed well with observations in the dam area of the TGR. Most of the deposits were found at the bottom of the main channel, whereas only a few deposits remained on the bank slopes. Due to the gravity-driven ﬂow of ﬂuid mud, the deposits that gathered in the deep channel formed a nearly horizontal surface. By considering the gravity-driven ﬂow, the averaged error of deposition thickness along the thalweg decreased from − 13.9 to 2.2 m. This study improved our understanding of the mechanisms of ﬁne-grained sediment transport in large reservoirs and can be used to optimize dam operations. a new method was developed and used to simulate the sedimentation pattern in the dam area of the TGR. The results showed that sediment would be deposited relatively uniformly in the dam area without considering the gravity-driven ﬂow. When the gravity-driven ﬂuid ﬂow of the ﬁne-grained deposits is considered, the distribution of the sediment deposits was more consistent with the observations than the results that do not consider gravity-driven ﬂow. The simulation results showed that the deposits mainly accumulated in the deep channel with a ﬂat surface, which agreed well with the observations. This study revealed the key role of the gravity-driven ﬂow of deposits in simulating the deposition of ﬁne-grained sediment in large reservoirs.


Introduction
The impoundment of large reservoirs results in numerous changes, including decreases in water flow velocity, sediment deposition, and sorting [1]. In the reservoir, coarser sediments tend to be deposited in the upstream regions of the reservoir, advancing downstream steadily but slowly in the form of a delta, whereas fine sediments reach the reservoir in the form of suspended load, by means of homogeneous or stratified dense flow [2].
Sedimentation is an important issue that receives a significant amount of attention in reservoir construction and operations [3]. The excess of sedimentation in a reservoir leads to sediment entrainment in waterway systems and hydraulic schemes [4]. During the demonstration phase and early operational periods of the Three Gorges Reservoir (TGR), several hydraulic research institutes studied the problem of reservoir sedimentation [5,6]. After the impoundment, the sediment influx conditions of the TGR changed significantly, and the actual reservoir sedimentation conditions were different from the predictions made during the demonstration phase [3,7]. With the increase in the water level in front of the dam during flood seasons, the trap efficiency increases [8,9]. In the TGR, deposited sediment is mostly fine sediment with median diameters of less than 0.01 mm, which was not forecasted by previous studies [7]. In addition, fine-grained sediments with median Considering the sedimentation features in the dam area of the TGR, we proposed a physical-based numerical model for fine-grained sediment deposition in large reservoirs. A method of modeling fine-grained sediment movement that includes gravity-driven motion in reservoirs was developed and applied, along with a 3D Reynolds-averaged numerical model of water flow and sediment transport, to simulate the uncommon sedimentation pattern in the dam area of the TGR. Our method can be applied to other reservoirs with similar sedimentation problems.

Features of Fine-Grained Sediment Deposition in Deep Reservoirs
The Three Gorges Dam (TGD), which is located in the Three Gorges reach along the trunk of the Yangtze River in China, is a key backbone project for flood control, navigation, and water resources development and has a total installed capacity of 18.2 GW and an average annual output of 84.7 billion kWh. The dam is located in Sandouping, 40 km upstream from Yichang City, as shown in Figure 1. The normal pool level of the Three Gorges Reservoir (TGR) is 175 m (relative to the Chinese water level zero point), which corresponds to an aggregated storage capacity of 39.3 billion m 3 . The TGR has been in operation since June 2003. The construction of large hydraulic projects such as the TGD has generated problems in navigation, environment, and ecology, all of which are related to sediment transport processes, including advection, diffusion, deposition, and erosion. Observations of the sedimentation pattern in the dam area of a typical section (S34#, see Figure 2 for the location) of the TGR have shown that the sediment mainly accumulated in the deep channel with a nearly horizontal surface (Figure 3) in the early years of reservoir operation. Conversely, deposits of fine-grained sediment in reservoirs typically have a relatively uniform distribution across the section. The physical and dynamic characteristics of the sediment deposits in the dam area can be used to understand this difference. Observations of the sedimentation pattern in the dam area of a typical section (S34#, see Figure 2 for the location) of the TGR have shown that the sediment mainly accumulated in the deep channel with a nearly horizontal surface (Figure 3) in the early years of reservoir operation. Conversely, deposits of fine-grained sediment in reservoirs typically have a relatively uniform distribution across the section. The physical and dynamic characteristics of the sediment deposits in the dam area can be used to understand this difference.
Before impoundment, the median size of the riverbed materials near the dam was approximately 0.10 mm. With the initial impoundment of the TGR, the water level rose to 135 m, and large amounts of siltation occurred in the river reach before the dam. In October 2003, the observed median particle size of the bed materials had decreased significantly to 0.02 mm. With further deposition, the median size had further decreased to 0.008 mm in May 2006. After the 156 m water storage operation of the TGR, the observed median size of the bed materials in the dam area was 0.005 mm in December 2006 and October 2007 [10]. The silt particles deposited in the dam area are extremely fine and require a long time to compact. Before compaction, the low-density deposits (the bulk density is approximately 1.00-1.20 t/m 3 ) in the form of sludge are called fluid mud. The fluid mud has significant fluidity and can flow or deform under its own weight. If the surface slope exceeds a certain value, the fluid mud will move down the slope. This may be the main reason that most of the deposits were found in the deep channel of the TGR and that the surface of the sludge is horizontal.  Before impoundment, the median size of the riverbed materials near the dam was approximately 0.10 mm. With the initial impoundment of the TGR, the water level rose to 135 m, and large amounts of siltation occurred in the river reach before the dam. In October 2003, the observed median particle size of the bed materials had decreased significantly to 0.02 mm. With further deposition, the median size had further decreased to 0.008 mm in May 2006. After the 156 m water storage operation of the TGR, the observed median size of the bed materials in the dam area was 0.005 mm in December 2006 and October 2007 [10]. The silt particles deposited in the dam area are extremely fine and require a long time to compact. Before compaction, the low-density deposits (the bulk density is approximately 1.00-1.20 t/m 3 ) in the form of sludge are called fluid mud. The fluid mud has significant fluidity and can flow or deform under its own weight. If the surface slope exceeds a certain value, the fluid mud will move down the slope. This may be the main   Before impoundment, the median size of the riverbed materials near the dam was approximately 0.10 mm. With the initial impoundment of the TGR, the water level rose to 135 m, and large amounts of siltation occurred in the river reach before the dam. In October 2003, the observed median particle size of the bed materials had decreased significantly to 0.02 mm. With further deposition, the median size had further decreased to 0.008 mm in May 2006. After the 156 m water storage operation of the TGR, the observed median size of the bed materials in the dam area was 0.005 mm in December 2006 and October 2007 [10]. The silt particles deposited in the dam area are extremely fine and require a long time to compact. Before compaction, the low-density deposits (the bulk density is approximately 1.00-1.20 t/m 3 ) in the form of sludge are called fluid mud. The fluid mud has significant fluidity and can flow or deform under its own weight. If the surface slope exceeds a certain value, the fluid mud will move down the slope. This may be the main  Therefore, the sedimentation process of fine sediment in deep areas of reservoirs can be summarized as several physical processes, including the dynamic process of sediment transport and settlement in flowing water, gravity-driven fluid flow of deposits, and compaction and consolidation.
The gravity-driven fluid flow of fine deposits is similar to that of fluid mud. Suppose that the surface slope of the fluid mud layer is J c . If J c > τ B /(γ s ∆z) (where τ B is the Bingham ultimate shear force, and γ s and ∆z are the weight and thickness of the fluid mud layer, respectively), the fluid mud would move down the slope under its own weight, as shown in Figure 4. Field observations have shown that the deposits in the deep channel are thick and that their surface is nearly horizontal in the dam area of the TGR (as shown in Figure 3). The main reason for this phenomenon is that the suspended sediment that was deposited during the early period of the TGR impoundment was so fine (median size of 0.004-0.01 mm) that the resulting deposits had a low initial dry bulk density and could move like fluid mud. The mud-like deposits moved downward of the surface slope due to gravity and accumulated in the deep channel of the riverbed, forming an uncommon sedimentation pattern in the early stage of the impoundment of the large reservoir.
The gravity-driven fluid flow of fine deposits is similar to that of fluid mud. Suppose that the surface slope of the fluid mud layer is . If / (where is the Bingham ultimate shear force, and and are the weight and thickness of the fluid mud layer, respectively), the fluid mud would move down the slope under its own weight, as shown in Figure 4. Field observations have shown that the deposits in the deep channel are thick and that their surface is nearly horizontal in the dam area of the TGR (as shown in Figure 3). The main reason for this phenomenon is that the suspended sediment that was deposited during the early period of the TGR impoundment was so fine (median size of 0.004-0.01 mm) that the resulting deposits had a low initial dry bulk density and could move like fluid mud. The mud-like deposits moved downward of the surface slope due to gravity and accumulated in the deep channel of the riverbed, forming an uncommon sedimentation pattern in the early stage of the impoundment of the large reservoir.

Numerical Simulation Method for Transport and Settlement of Fine-Grained Sediment
Based on the analysis presented above, the whole processes of fine-grained sediment in the deep reservoir can be summarized as transport and settlement, gravity-driven fluid flow, and consolidation. The 3D Reynolds-averaged numerical model that depicts flow, sediment transport, and sediment settlement in a reservoir can be found in the related literature [20]. Here, we briefly describe the method.
The governing equations for water flow are given as:

Numerical Simulation Method for Transport and Settlement of Fine-Grained Sediment
Based on the analysis presented above, the whole processes of fine-grained sediment in the deep reservoir can be summarized as transport and settlement, gravity-driven fluid flow, and consolidation. The 3D Reynolds-averaged numerical model that depicts flow, sediment transport, and sediment settlement in a reservoir can be found in the related literature [20]. Here, we briefly describe the method.
The governing equations for water flow are given as: where u i (i = 1,2,3) are velocity components; f i is the gravitational force per unit volume; ρ is the density of water; and p is the pressure. In this model, turbulence stresses τ ij are calculated with the standard k − ε turbulence model. The 3D nonuniform suspended sediment transport equation can be expressed in Cartesian coordinates as: ∂s k ∂t where u, v, and w are the flow velocities in the x, y, and z directions, respectively; s k is the sediment concentration of group k; ω sk is the settling velocity of group k; and ε s is the sediment diffusion coefficient. The equation of bed evolution caused by suspended sediment transport (deposition) is where γ s is the dry bulk density of the deposits; s kb and s kb * are the real and equilibrium sediment concentrations of group k near the bed, respectively; Z b is the bed level; and N s is the total number of sediment groups.

Numerical Simulation Method for Gravity-Driven Motion of Fine-Grained Sediment Deposits
Solid consolidation can be determined using existing methods. Here, we mainly discuss the method of simulating the gravity-driven flow of the fine-grained sediment deposits. Newly deposited fine-grained sediment deposits are a non-Newtonian fluid and are generally considered to be a Bingham fluid. Before movement, it is necessary to overcome the Bingham ultimate shear force. In this paper, a critical slope is used as a criterion for the destabilizing flow of fine-grained sediment deposits. The shallow water flow model is used to describe the gravity-driven flow process of the destabilized deposits.
When the surface slope of fine-grained sediment deposits J c satisfies the condition J c > τ B /(γ s ∆z) (where τ B is the Bingham ultimate shear force, and γ s and ∆z are the weight and thickness of the fluid mud layer, respectively), the deposits will move down the slope under their own weight [24]. Fine-grained deposits with different bulk densities have different Bingham ultimate shear forces. Here, we apply the formula for the Bingham ultimate shear force on fluid mud on a muddy shore [25], where α and β are empirical parameters; α = 9 × 10 −9 ; and β = 16.719 based on previous research [25]. After losing stability, the fine-grained deposit will flow. As the thickness of the deposits is small compared to the horizontal scale of the reservoir, we assume that the flow of the fine-grained deposit can be described by the 2D shallow flow equations.
Continuity equation: ∂Z ∂t Momentum equation: ∂N ∂t where Z is the surface elevation of the fine-grained deposits; M = uh; N = vh; h is the thickness of the deposits; u and v are the velocities in the x and y directions, respectively; q is the source/sink of the fine deposit in a unit area, which can be calculated according to the instability as q = (∆z − ∆z c )/dt, where ∆z c is the maximum thickness that can remain stable, and ∆z c = τ B (J c γ s ) ; g is the acceleration due to gravity; D is the kinematic viscosity; and n is the integrated Manning's roughness. As gravity is the main driving force of the flow, the diffusion terms in Equations (7) and (8) can be neglected.

Numerical Simulation Procedure
Consistent with the 3D flow and sediment transport model, the gravity-driven flow equations of fine deposits are discretized based on the finite volume method, which has the advantage of conservation. The SIMPLEC algorithm is used to solve the problem. The solution process is shown in Figure 5 and illustrated as follows.
sediment deposition and to obtain the distribution of the deposition thickness .
Step (2): Based on the silt thickness and terrain slope, determine whether the deposits in a cell are stable or not. If / , the deposit is unstable and the source term q will be calculated; otherwise, the deposit is stable and q is zero.
Step (3): Solve the deposit flow equations and apply the SIMPLEC algorithm to obtain the surface elevation of the deposits after being redistributed.

Results
Based on the topographic, hydrological, and sediment data measured in the initial stage of the TGR impoundment, a 3D Reynolds-averaged water flow and sediment transport numerical model in the dam area of the TGR (from Miaohe to the dam, see Step (1): Use the 3D numerical model to simulate the flow, sediment transport, and sediment deposition and to obtain the distribution of the deposition thickness ∆z.
Step (2): Based on the silt thickness and terrain slope, determine whether the deposits in a cell are stable or not. If J c > τ B /(γ s ∆z), the deposit is unstable and the source term q will be calculated; otherwise, the deposit is stable and q is zero.
Step (3): Solve the deposit flow equations and apply the SIMPLEC algorithm to obtain the surface elevation of the deposits after being redistributed.

Results
Based on the topographic, hydrological, and sediment data measured in the initial stage of the TGR impoundment, a 3D Reynolds-averaged water flow and sediment transport numerical model in the dam area of the TGR (from Miaohe to the dam, see Figure 2) was established. The simulated area was approximately 15 km long. The initial terrain for the simulation was measured in June 2003 ( Figure 6). Flow and sediment data at Miaohe from June 2003 to November 2008 were used for the inflow conditions. The water level in front of the dam was used for the outflow conditions. The numbers of grid points in the longitudinal, lateral, and vertical directions were 198, 97, and 18, respectively.     Figure 8 shows the calculated suspended load concentrations on typical cross-sections. The sediment concentration decreases much more rapidly in wider valley areas. Due to the impact of the secondary currents shown in Figure 7, the maximum concentration occurs at the convex side, whereas the minimum occurs at the concave side. Figure 9 shows a comparison of the observed and simulated velocities and suspended load concentrations on typical cross-sections (see Figure 2 for the locations of the sections) for a flow discharge of 24,600 m 3 /s. In general, the 3D numerical model accurately predicted the flow patterns and suspended load concentrations in the dam area of the TGR.

Features of Sedimentation Pattern
The river reach in front of the TGD is a wide valley of crystalline rocky hills with significant slopes on both banks. Fine-grained sediment settled in this area after the impoundment of the TGR. If the thickness of the deposits in a certain location exceeds a critical value, the deposits will flow down the slope under their own weight. Based on the verifications of the flow velocity and sediment concentration (Figure 9), Figure 10 shows the sediment thickness near the dam calculated by traditional 3D flow and sediment transport modeling (without considering the gravity-driven fluid flow of the deposits), and Figure 11 shows the results calculated using the proposed method. Although the total amounts of sediments are consistent, the sediment distributions are quite different. The sediment thickness calculated by the traditional method, which does not consider the fluid flow of the fine-grained deposit, is relatively uniform in map view ( Figure 10). In contrast, the proposed model ( Figure 11) predicts that most of the deposits will accumulate in the deep channel and that a few deposits will remain on the flat terraces on both sides, which is consistent with the observations. Figure 12 shows comparisons between the observed and calculated sedimentation on several typical cross-sections (section locations are shown in Figure 5). The cross-sectional sedimentation pattern that considers gravity-driven motion is more similar to the observations than the pattern that does not consider the motion. This again confirms the significance of the gravity-driven flow of fine deposits.

Features of Sedimentation Pattern
The river reach in front of the TGD is a wide valley of crystalline rocky hills with significant slopes on both banks. Fine-grained sediment settled in this area after the impoundment of the TGR. If the thickness of the deposits in a certain location exceeds a critical value, the deposits will flow down the slope under their own weight. Based on the verifications of the flow velocity and sediment concentration (Figure 9), Figure 10 shows the sediment thickness near the dam calculated by traditional 3D flow and sediment transport modeling (without considering the gravity-driven fluid flow of the deposits), and Figure 11 shows the results calculated using the proposed method. Although the total amounts of sediments are consistent, the sediment distributions are quite different. The sediment thickness calculated by the traditional method, which does not consider the fluid flow of the fine-grained deposit, is relatively uniform in map view ( Figure 10). In contrast, the proposed model ( Figure 11) predicts that most of the deposits will accumulate in the deep channel and that a few deposits will remain on the flat terraces on both sides, which is consistent with the observations. 1, 13, x FOR PEER REVIEW 11 of 17 Figure 10. Computed sediment distribution using the traditional method.   . Computed sediment distribution considering gravity-driven processes. Figure 12 shows comparisons between the observed and calculated sedimentation on several typical cross-sections (section locations are shown in Figure 2). The crosssectional sedimentation pattern that considers gravity-driven motion is more similar to the observations than the pattern that does not consider the motion. This again confirms the significance of the gravity-driven flow of fine deposits.  Figure 13 compares the surface elevations along the thalweg between the observations and the results calculated with and without considering the gravity-driven deposit flow. It can be seen that the calculated results that consider the deposit flow match the observations quite well. Except the point at the entrance, the maximum error at the other points is 5.1 m, which is located at section S31# (see Figure 2 for the location). Based on measured data, the averaged and maximum deposition thickness along the thalweg are 22.

Difference between Calculation and Measurement
Although the comparison strongly confirmed the necessity and advantage of considering the gravity-driven movement of the deposits, there was some discrepancy between the observations and simulating results.
Near the dam, as shown in Figures 12a and 13, the depth of the deposit at the thalweg is greater than that of the observation. This is probably because this section is near the discharge hole, so the fluid deposits may be scoured and flushed by the flow. The effect of discharge flow was not simulated precisely in this study, which may have caused the deviation. Figure 14 shows the computed velocity distribution in front of the dam. In this work, the outlets of the power plant were generalized as a single outlet area, which would flatten the flow speed nearby. Another possible reason is the terrain flattening caused by mesh interpolation in the numerical simulation. Without terrain flattening, some of the deposit may be blocked by local terrain undulation, reducing the amount of deposits that gathered in the deep channel. It is easy to imagine that those cross-sections with wider sides, more undulating, relatively gentle side slopes would suffer greater impact from terrain flattening. As the reach of the near-dam area gradually widens from Miaohe to the dam, the terrain flattening error is generally more obvious at the sections before the dam, which can explain the relatively larger discrepancy within 5 km before the dam than that of the upstream results in Figure 13.

Difference between Calculation and Measurement
Although the comparison strongly confirmed the necessity and advantage of considering the gravity-driven movement of the deposits, there was some discrepancy between the observations and simulating results.
Near the dam, as shown in Figures 12a and 13, the depth of the deposit at the thalweg is greater than that of the observation. This is probably because this section is near the discharge hole, so the fluid deposits may be scoured and flushed by the flow. The effect of discharge flow was not simulated precisely in this study, which may have caused the deviation. Figure 14 shows the computed velocity distribution in front of the dam. In this work, the outlets of the power plant were generalized as a single outlet area, which would flatten the flow speed nearby. Another possible reason is the terrain flattening caused by mesh interpolation in the numerical simulation. Without terrain flattening, some of the deposit may be blocked by local terrain undulation, reducing the amount of deposits that gathered in the deep channel. It is easy to imagine that those cross-sections with wider sides, more undulating, relatively gentle side slopes would suffer greater impact from terrain flattening. As the reach of the near-dam area gradually widens from Miaohe to the dam, the terrain flattening error is generally more obvious at the sections before the dam, which can explain the relatively larger discrepancy within 5 km before the dam than that of the upstream results in Figure 13.

Shortcomings Of Present Method
In this work, shallow water equations are applicable to a Newtonian fluid. For a Bingham fluid, when integrating the continuity and momentum equations with depth to obtain the 2D governing equations, the nonlinear relationship between the shear stress and strain would generate significant differences with the shallow water equations [26]. If the shallow water equations are used, the resistance term would be more complicated than that of Manning's formula [27]. Therefore, the modeling of gravity-driven flow of finegrained deposits should be improved further. In addition, due to the discretization, the topographic data used in the numerical simulation are not consistent with the real topography. This also contributed to the discrepancy between the simulated results and observations. However, considering the gravity-driven fluid flow of the deposits, the simulation results of the overall sedimentation trend and distribution showed good agreements with the observations.

Conclusions
The deposition processes of fine-grained sediment in large reservoirs can be summarized as the settling of suspended sediment, gravity-driven fluid flow of the deposits, and sludge consolidation and compaction. Based on the physical mechanisms, we analyzed the dynamic causes of the patterns of fine-grained sedimentation in the dam area of the TGR during the initial stage of impoundment. The sediment deposited in the dam area of the TGR during this stage is mostly fine-grained sediment with a grain size generally between 0.004 and 0.010 mm. These fine-grained sediment deposits have a low dry bulk density and can move downslope toward the bottom of the riverbed due to the effect of gravity. These fluid deposits finally gather in the deep channel and have a horizontal surface, which is an uncommon feature of sedimentation in large reservoirs. The rapid increase in the surface elevation of the silt in front of the dam threatens the operations of the dam.

Shortcomings of Present Method
In this work, shallow water equations are applicable to a Newtonian fluid. For a Bingham fluid, when integrating the continuity and momentum equations with depth to obtain the 2D governing equations, the nonlinear relationship between the shear stress and strain would generate significant differences with the shallow water equations [26]. If the shallow water equations are used, the resistance term would be more complicated than that of Manning's formula [27]. Therefore, the modeling of gravity-driven flow of fine-grained deposits should be improved further. In addition, due to the discretization, the topographic data used in the numerical simulation are not consistent with the real topography. This also contributed to the discrepancy between the simulated results and observations. However, considering the gravity-driven fluid flow of the deposits, the simulation results of the overall sedimentation trend and distribution showed good agreements with the observations.

Conclusions
The deposition processes of fine-grained sediment in large reservoirs can be summarized as the settling of suspended sediment, gravity-driven fluid flow of the deposits, and sludge consolidation and compaction. Based on the physical mechanisms, we analyzed the dynamic causes of the patterns of fine-grained sedimentation in the dam area of the TGR during the initial stage of impoundment. The sediment deposited in the dam area of the TGR during this stage is mostly fine-grained sediment with a grain size generally between 0.004 and 0.010 mm. These fine-grained sediment deposits have a low dry bulk density and can move downslope toward the bottom of the riverbed due to the effect of gravity. These fluid deposits finally gather in the deep channel and have a horizontal surface, which is an uncommon feature of sedimentation in large reservoirs. The rapid increase in the surface elevation of the silt in front of the dam threatens the operations of the dam.
By coupling a traditional 3D water flow and sediment transport model with a newly proposed numerical method for the gravity-driven flow of fine-grained sediment deposits, a new method was developed and used to simulate the sedimentation pattern in the dam area of the TGR. The results showed that sediment would be deposited relatively uniformly in the dam area without considering the gravity-driven flow. When the gravity-driven fluid flow of the fine-grained deposits is considered, the distribution of the sediment deposits was more consistent with the observations than the results that do not consider gravity-driven flow. The simulation results showed that the deposits mainly accumulated in the deep channel with a flat surface, which agreed well with the observations. This study revealed the key role of the gravity-driven flow of deposits in simulating the deposition of fine-grained sediment in large reservoirs.