Abstract
To evaluate the evolutionary processes guiding the formation of the tailings-water mixtures produced by the instantaneous collapse of tailings ponds and the influence of these on downstream facilities, a 2D simulation model with reasonable boundary and working conditions derived from actual engineering practice was built in this study, and the relationship between dam-break elevation and impact on downstream facilities was also analyzed to determine the relevant mechanism of influence. Computational results indicated that lowering the dam-break elevation caused the maximum velocity and flooding depth, along with the flooded area at monitoring points, to gradually increase. The occurrence times of maximum velocity and flooding depth were also gradually moved forward as the breaking elevation was reduced; this effect is directly related to the increase in the total potential energy at the lower break elevations. Further simulations of sand-prevent dams with different heights located downstream from a tailings pond were carried out to identify methods for mitigating the impact of dam failure. The results revealed that increasing the height of the sand-prevent dam reduced the production of tailings mixtures. Based on the results, the construction of a sand-prevent dam with a crest elevation equal to that of the starter dam was recommended.
1. Introduction
Tailings ponds are specialized structures for the deposition of tailings produced by mining extraction [1]. For historical and policy purposes, tailings ponds in China can be characterized as being numerous but with small capacities and have often been constructed with little attention to design, with dam structures that do not meet stability requirements. The water levels in the reservoirs of upstream tailings ponds are prone to change when there is heavy rain or impervious structural failure [2–4]. Because of the complex physical and mechanical properties of tailings and their high levels of physical discreteness, the risk of dam break is extremely high when catastrophic floods [5] or earthquakes occur. Dam breaks in tailings ponds can result in sudden and significant damage [6], especially when there are facilities or densely populated villages downstream, and can cause serious consequences—such as casualties, property loss, and environmental pollution—that are difficult to remedy.
As the causes of tailings pond failure are numerous and complex, researchers have had to systematically investigate and analyze them. For example, Mei and Wang [7] conducted the statistical analysis and countermeasure research on the causes of tailings dam failure accidents in China. Berghe et al. [8] systematically analyzed the risks of tailings dam breakage; based on a Chinese case study, Wei et al. [9] proposed specific requirements for the design, construction, and management of tailings ponds. Dobry and Alvarez [10] analyzed post-earthquake dam-break accidents and found that failure typically occurred during or immediately after the earthquake and was induced by slope slide and flow. Robinson and Toland [11] discussed case histories covering a wide variety of embankment and foundation conditions and noted that seepage problems associated with foundation conditions have a pronounced impact on the overall ability of facilities to retain their tailings perimeters. Shakesby et al. [12, 13] studied the failure and subsequent development of gold mine tailings dam flow slides based on a combination of multiple factors, including poor basal drainage, too-steep perimeter walls, wall saturation, and the effects on tailings saturation of the accumulation of basal sediments through the continued spitting of slurry during periods of heavy rainfall. Priulli et al. [14] presented a refined computing model for simulating rapid-flow movement across a 3D terrain. Rico et al. [15] carried out a detailed search and reevaluation of known historical cases of tailings dam failure in Europe and found that failure was most often related to unusual rainfall and that more than 90% of incidents occurred inactive mines, with only 10% occurring in abandoned ponds. Their results indicated an urgent need for EU regulations covering the technical standards of tailings disposal. Fourie et al. [16] assessed the 1994 Merriespruit, Virginia gold tailings dam failure, which resulted in 17 deaths, and found that it was produced by a large volume of metastable tailings and the overtopping and erosion of the impoundment wall, which combined to produce static liquefaction of the tailings and a consequent flow failure. Zhang et al. [17, 18] carried out an experimental model study on dam breakage and the tailings pond evolution law and proposed a model test method for predicting and preventing the dam-break process. Li [19] analyzed the influence of rainfall factors and safety warning technology on tailings dam breakage and, based on physical and numerical simulation of a dam example, determined the evolution rule and water flow characteristics for instantaneous and progressive tailings dam break. Mei [20] researched the mechanism of and online early-warning methods for tailing dam failure, and Xu [21] studied the dynamic process and post-failure effect of the overtopping failure of a tailings dam.
To evaluate the risk level of tailing dam failure accurately, Wang et al. [22] established a matter-element extension evaluation model. Considering that the stability of a tailings dam is affected by many random and fuzzy factors that are difficult to quantitatively estimate, Li et al. [23] conducted fuzzy theory-based research on tailings dam failure risk evaluation models and proposed a fuzzy theoretical model for the risk assessment of tailings ponds. Similarly, Liang et al. [24] proposed a dam-break risk assessment model based on variable weight synthesis and analytic hierarchy processes. Chen et al. [25] proposed a dam-break risk classification model based on interpretive structural modeling and factor frequency methods. Xiong [26] established the application of downstream dam-break risk fields to tailings ponds based on the use of geographic information systems and risk control theory. Zhang [27] also developed a risk assessment model and control strategy for tailings dam failure. Model testing and numerical modeling are both important tools for evaluating the serious consequences of tailings pond dam failure, and researchers including Li et al. [28], Li et al. [29], Wei et al. [30], Li [31], Liu [32], Lu [33], and Qin [34], have simulated tailings pond dam breaks under different working conditions. Nevertheless, there remains a lack of research on the impact of dam breakage on downstream areas.
If a dam is broken for any reason, its accumulated tailings can escape and flow downstream along the existing terrain. In this study, the evolution law of the spreading range, flow rate, and submergence depth of the mixture of tailings and water under different dam-break elevations was determined using a 2D numerical simulation model of dam failure based on the Mike 21 hydrodynamics software package. A practical terrain- and dam data-based model of a tailings dam was established and, using a set of reasonable river bed resistance coefficients, an evaluation of the effects of tailings pond release on downstream facilities was carried out.
2. Governing Equations of the Dam-Break Simulation
The sand flow caused by tailings discharge after dam collapse is inherently a landslide or debris flow in nature and can be assumed to be a special form of movement between “fluid” and “granular body.” This flow can be described by the dynamic equation and continuity equation similar to the fluid flow.
Therefore, this study is based on the following assumptions. The tailings sand is an isotropic continuum; the flow of the tailings sand conforms to the Bingham flow mode; the relationship between shear stress and shear strain is linear but does not pass through the origin; at a certain yield stress, the collapse of the tailings dam occurs under the condition that the dam is piled high. Thus, the water storage in the tailings dam corresponding to the overtopping water level of the tailings dam is the largest.
Based on a solution of the 2D incompressible average Navier–Stokes equation established using the Boussinesq and two-dimensional incompressible assumptions, the sand and water escape from the breakdown of a tailings pond can be approximately described by the following 2D average depth equations [35]: Continuity equation: Momentum equations:where t is time, x and y are Cartesian spatial coordinates, is the height of the water, d is still water depth, is the total water depth, is the Coriolis parameter (where is the angular velocity of rotation, is the latitude), is the acceleration of gravity, is the density of water, , , , and are the respective components of radiation, is the reference density of water, S is the flow size of the point source, and and are the velocity components in the x- and y-directions, respectively, which are determined to the average velocity by water depth as follows:
The lateral stress component comprehensively considers the influence of viscous and turbulent friction, differential advection, etc., whereas the vortex viscosity formula is estimated based on the average velocity gradient:
3. Dam-Break Simulation
3.1. Basic Information on the Tailings Pond Model
The starter dam modeled for this study was a clay core dam with bottom and crest elevations of 1,156.0 and 1,180.0 m, respectively, and an external slope ratio of 1 : 2.0. The final accumulation elevation of the dam crest was 1,195.0 m and the external slope ratio of the accumulation dam was 1 : 3.0. The tailings pond was a fourth-class pond with a total capacity of approximately 629,000 m3. Several facilities such as villages and high-speed railway lines were located approximately 510 m downstream of the tailings pond to the southeast (see Figure 1 for details). The main cross section of the tailings pond is shown in Figure 2.
A flow chart of the tailings pond dam-breaking simulation is shown in Figure 3. The scope of calculation was determined based on the collection of a sufficient amount of survey, topography, and other basic data. The reliability, consistency, and representativeness of the data were then evaluated. The numerical scheme used to model the dam-break simulation was determined based on an analysis of reasons for dam break and outburst shape, duration, and elevation. Using the resulting calculation scheme, the terrain was digitized, the grid was divided to meet the requirements of computational accuracy, and reasonable calculation time series and boundary conditions were set. Finally, the resulting flow rate, inundation range, and depth of discharge were analyzed to obtain the effects of dam break on the downstream facilities.
3.2. Grid Division, Boundary Conditions, and Location of the Monitoring Points
A finite volume method based on the mesh center was used to spatially discretize the governing equations. The horizontal plane was covered in an unstructured mesh comprising a mixed grid of triangles and quadrilaterals to provide an optimal fitting to the complex geometric topography through an arrangement of small and large elements in key and non-key areas, respectively, while smoothing the boundary. The grid subdivision of the dam-break simulation, which is shown in Figure 4, adopted a grid size of 50 m and divided the study area into 59,061 elements. As all of the boundaries were closed, there was a zero flow velocity in the vertical direction.
To analyze the flow velocity and deposition of discharged material at the inclined shaft and drain opening during the discharge process, monitoring points were set at these locations to observe the changes in discharged material over time. The locations of the monitoring points are shown in Figure 5.
The break patterns of the tailings dam can be generally divided into total and local failure in terms of scale and instantaneous and gradual dam failure in terms of time. As the impact force of the flow causing the dam break is strong, the entire process from the beginning of the break to the formation of a stable rupture section is very short; therefore, the dam break can be treated as instantaneous and total collapse for safety reasons.
The amount of tailings and water discharged is directly related to the amount of water stored in the tailings pond. In this study, the storage capacity characteristic, corresponding to the accumulation elevation of +1195 m of the tailings dam, was set as the initial condition of the simulation, whose water storage amount was the largest. Thus, the debris flow formed by the dam break was also the largest. The crossing line between the tailings pond and mountain was defined as a dike, whose crest elevation was the current elevation. The tailings dam was simulated by the axis of the highest accumulation dam and defined as a dike with a constant elevation of the crest, whose elevation adopted the dam-break elevation. The sand-prevent dam used in this study was also defined as a dike, whose crest elevation was also set as a constant value. The other boundary conditions adopted current topographic data. A total of 60 time steps was set in this simulation, and each of the time step was 30 s.
3.3. Bed Resistance Parameter
A variation of Manning’s coefficient was used as a shear force function to simulate the bed resistance parameter in modeling the different flow characteristics of tailings slurry and water. The relationship between the shear force and Manning’s coefficient is shown in Figure 6.
In Figure 6, τc is the critical shear force as the tailings begin to flow, τcu is the critical shear force of the tailing slurry as it transitions from a yield pseudoplastic body to a Bingham plastic body, and Ml is the minimum Manning coefficient, which represents the bottom friction as the slurry begins to flow; here, M = 1/n, where n is the roughness coefficient. There are several facilities in the downstream riverbed and the mountains on both sides are steep, with a few trees on the riverbed. The roughness values in Table 1 were used to select a minimum Manning coefficient of 3.2. A maximum Manning coefficient, which represents the bottom friction as the slurry transitions from a pseudoplastic body to a Bingham plastic body, of Mu = 32, was selected.
To assess the impact of dam breakage elevation on the downstream facilities, three different breaking elevations were modeled—two-thirds of the accumulation dam height, one-third of the accumulation dam height, and the elevation of the dam crest of the starter dam—and the resulting discharge evolution processes and distributions of tailings deposition were compared and analyzed.
3.4. Analysis of Tailings Evolution Results
3.4.1. Break Occurring Immediately above the Two-Thirds of Accumulation Dam Height Mark
Because an instantaneous dam break occurs quickly and the maximum flow appears in the initial stage of the break, the scope of spread within 30 min of breaking was analyzed. The corresponding elevation at two-thirds of the dam height was +1,190 m. The evolution of the tailings-and-water mixture flow following an instantaneous tailings pond collapse above this altitude is shown in Figure 7. In this case, the flow spreads downstream along the valley under the action of accumulated potential energy. 8.5 minutes after the dam break, the mixture of tailings and water reaches the inclined shaft opening and continues to spread downstream; 9.5 min after the dam break, tailings, and water flood the drain opening; 30 min after the break, the reservoir has ceased discharging. Although the quantity of mixed tailings and water at the inclined shaft is small, it is large at the drain opening; this effect is directly related to the respective elevations and topographic conditions of the two locations.
(a)
(b)
(c)
(d)
(e)
(f)
3.4.2. Instantaneous Breakage of the Dam at One-Third of Dam Height
The evolution process of tailings and water flow following a break at one-third of the dam height (elevation +1,185.0 m) is shown in Figure 8. The working conditions are similar to those occurring following the break at +1,190.0 m: 4.5 min after the dam break, the mixture of tailings and water reaches the inclined shaft; 6.5 min after the break, it reaches the drain opening; 30 min after the break, the discharge stops, and a large number of tailings and water have accumulated at the drain opening. The quantity of mixture accumulating at the inclined shaft, however, is less than that under the two-thirds case and does not affect its primary functionality.
(a)
(b)
(c)
(d)
(e)
(f)
3.4.3. Instantaneous Breakage at Starter Dam Crest Height
Finally, the breakage of the dam at the elevation of the crest of the starter dam was modeled. The evolution process of tailings-and-water flow is shown in Figure 9. Four minutes after the dam breaks, the mixture reaches the inclined shaft; 5.5 min later, the mixture reaches the drain opening.
(a)
(b)
(c)
(d)
(e)
(f)
3.5. Analysis of Flooding Depth and Flow Velocity
Under a basic assumption of our dam-break simulation, the total amount of discharge will be equivalent to the reservoir capacity above the elevation of the break. The changes in flooding depth and flow velocity, respectively, overtime at the monitoring points in front of the inclined shaft and drain opening downstream of the tailings pond are shown in Figures 10 and 11.
(a)
(b)
(a)
(b)
It is seen that the flooding depths at both monitoring points first increase to maximum values and then decrease gradually over time. As the dam-break elevation is reduced, the maximum flooding depth increases, and the maximum value is reached within a shorter time. This is directly related to the water-capacity curve of the tailings pond (i.e., the shape of the valley in which the pond is located). Reducing the dam-break position increases the total potential energy of the discharged material, and therefore the total energy of the discharged material, enhancing both the velocity and volume of flow. The maximum flooding depth at the drain opening exceeds that at the inclined shaft because the inclined shaft is located at a higher elevation than the drain opening.
The flow velocity curves in front of the inclined shaft and drain opening follow similar laws in terms of the forms of their respective flooding depth curves: as the elevation of dam break decreases, the maximum velocity gradually increases, and the time at which the maximum occurs gradually advances. After reaching its maximum, the velocity decreases gradually over time and finally settles into its initial state. In the second halves of the velocity curves, corresponding to the drain opening, the fluctuation in velocity decreases owing to the proximity of the closed boundary downstream of the drain opening. When the discharge reaches this boundary, it flows back under the action of boundary resistance, resulting in a decrease in the fluctuation of the velocity curve.
4. Influence of Sand-Prevent Dam on the Amount of Discharge
The tailings pond modeled in this study was an overhead reservoir with villages and high-speed railway facilities present within 1 km downstream. According to the analysis results described in Section 3, the dam-break tailings-and-water mixture would flood this infrastructure. To reduce the impact of discharge on the villages and high-speed railway facilities, we examined the effects of placing a sand-prevent dam, built using the same material used to construct the starter dam, 150.0 m downstream of the axis of the starter dam.
Three different sand-prevent dam schemes, with crest elevations of 1,180.0, 1,170.0, and 1,160.0 m respectively, were modeled. According to the analysis results presented in Section 3, the impact on downstream facilities would be most severe following a tailings dam break at the elevation of the starter dam crest. Therefore, for each sand-prevent model, we analyzed the evolution rules of discharge flow and submergence depth following a tailings dam break at an elevation of 1,180.0 m.
4.1. Analysis of Flooding Depth and Flow Velocity
The discharge flow in front of the inclined shaft is shown in Figure 12. It is seen from the figure that the discharge flow evolution law is essentially the same at all three sand-prevent dam crest elevations; in each case, there is a steep increase from the initial flow value of zero and a gradual decrease following the peak of the curve. As the height of the sand-prevent dam increases, the curve decreases more sharply following the peak. At a sand-prevent dam crest elevation of 1,180.0 m, the discharge flow velocity at the inclined shaft is only 8.0% of that occurring without a sand-prevent dam.
(a)
(b)
(c)
(d)
The flooding depths in front of the inclined shaft under the respective sand-prevent dam models are shown in Figure 13. It is seen that the submergence depth evolution rules under the different sand-prevent dam schemes are also essentially the same; following a steep increase from an initial depth of zero, the depth peaks and then gradually decreases. For the sand-prevent dam with a height of 1,180.0 m, the maximum submergence depth in front of the inclined shaft is 0.13 mm, which is only 1.0% of the depth occurring without a sand-prevent dam; this indicates that the sand-prevent dam can intercept most of the discharge under this working condition.
(a)
(b)
(c)
(d)
The flow velocities at the drain opening under the respective sand-prevent dam models are shown in Figure 14. It is seen that the evolution rules are essentially the same, with each showing a steep increase from an initial velocity of zero and a gradual decrease in fluctuation after the peak. Increasing the height of the sand-prevent dam causes the curve to decrease more sharply after the peak. At a sand-prevent dam crest elevation of 1,180.0 m, the flow velocity at the drain opening is only 22.1% of that obtained without a sand-prevent dam.
(a)
(b)
(c)
(d)
The flooding depths at the drain opening under the respective sand-prevention models are shown in Figure 15, from which it is evident that the flooding depths at the drain opening follow essentially the same evolution law; in each case, there is a steep increase from an initial depth of zero and then a gradual decrease after the peak. At a sand-prevent dam crest elevation of 1,180.0 m, the maximum flooding depth at the drainage outlet is only 38.9% of that obtained without a sand-prevent dam, indicating that the sand-prevent dam can intercept most of the discharged tailings and water under this working condition.
(a)
(b)
(c)
(d)
4.2. Analysis of Discharge Velocity and Flooding Depth under Different Sand-Prevent Dam Schemes
The maximum flooding depths at the inclined shaft and drain opening under each sand-prevention scheme are shown in Figure 16. As the height of the sand-prevent dam increases, the flooding depths at both positions gradually increase, indicating that the flooding depth near the downstream facilities will be directly related to the height of the sand-prevent dam. Increasing the height of the sand-prevent dam reduces the flooding depth, and, therefore the impact of a tailings dam break on the downstream area.
(a)
(b)
The maximum values of discharge velocity at the inclined shaft and drain opening under each sand-prevent dam scheme are shown in Figure 17. As the height of the sand-prevent dam increases, the discharge velocities at both locations gradually increase, indicating that the discharge velocities at the facilities downstream will be directly related to the height of the sand-prevent dam. Increasing the dam height will reduce the discharge flow velocity, thereby reducing the impact of a tailings pond dam break on the downstream facilities.
(a)
(b)
4.3. Discussion
According to the existing Chinese national codes or regulations for design or safety of tailings pond, there is no suggestion regarding the construction of sand-prevent dams to reduce the impact of dam-break accidents. Only tailings collection dams are recommended in centerline tailings ponds, which are typically constructed in the downstream direction of the tailings dam and used to block tailings entrained by rainwater erosion.
In China, a tailings dam with inhabited areas or facilities downstream within 1 km is defined as a high-risk tailings dam, referred to as the “on-head” tailings dam; therefore, 1 km is a controlling distance. The spread range of dam-break debris flow is closely related to the characteristics of the discharged materials, and the influence range is more than 1 km under extreme conditions. For example, on November 5, 2015, the Fundao tailings pond in Samarco Iron Mine in Brazil collapsed owing to a small earthquake, which triggered the liquefaction of the super-high dam that was already close to saturation. At least 19 people were killed, with an estimated 32 million m3 of tailings being discharged and flooding 158 houses in the village of Bento Rodrigues, 5 km downstream of the pond. On January 25, 2019, Brazil’s Feijiao tailings dam collapsed, releasing up to 9.8 million m3 of tailings, resulting in 270 deaths, with the main victims being residents of the Vila Ferteco community, 1 km away from the dam. The cases above show that, if no blocking measures are taken, dam collapse can have a significant impact on downstream areas, where residents or facilities exist. Moreover, from a safety perspective, distances beyond 1 km along the flow path should be studied further.
The theoretical research above shows that a sand-prevent dam located downstream of the tailings dam is an effective measure to reduce the consequences of dam collapse for residents or facilities that are downstream of the tailings pond. Therefore, from a safety perspective, a sand-prevent dam is suggested to be constructed to further ensure safety.
5. Conclusions
In this study, technical data obtained from an actual tailings dam project were used to quantitatively model the flooded area and flooding depth following dam breakage using numerical simulation. The impact of tailings pond dam failure on downstream facilities was analyzed, producing the following conclusions:(1)Five minutes after dam breakage, the tailings-water mixture flow will have flooded the downstream high-speed railway facilities and drain opening, of which the scope of influence is longer than 1 km. This will seriously influence the normal operation of the railway and damage downstream villages and facilities, and a portion of the discharge will pile up in the natural riverbed with detrimental ecological impact.(2)Numerical analysis of the dam-break process revealed that the evolution curves of flow velocity and flooding depth follow similar laws, with both increasing sharply, reaching a maximum value, and finally decreasing at a gradually slowing rate.(3)The evolution curves of flow rate and submerged depth at typical positions in front of the inclined shaft and drain opening indicate that reducing the dam-breaking position causes the maximum velocities and submerged depths at these positions to increase gradually. Both effects are related to the fact that the discharge capacity of the reservoir increases as the breakage height of the dam decreases. Furthermore, the times of occurrence of maximum velocity and submerged depth will gradually advance as the dam-break elevation decreases because reducing the height of breakage increases the total potential energy accumulated at the dam-break position.(4)Numerical modeling of the effects of placing sand-prevent projects with different heights revealed that increasing the height of the sand-prevent dam reduced the flooding depth and flow velocity of the discharge. For a sand-prevent dam with a crest height equal to that of the starter dam, the submerged depth and flow velocity of the tailings-and-water mixture in front of the inclined shaft will be only 1.0 and 8.0%, respectively, of the depth and velocity achieved without a sand-prevent dam. Constructing a sand-prevent dam of this height will reduce the flooding and flow velocity in front of the drainage outlet to only 38.9 and 22.1% of the respective values without a sand-prevent dam. These results indicate that a sand-prevent dam with a crest elevation similar to that of the starter dam should be set up downstream of the overhead tailings pond to enhance the safety of nearby facilities and significantly reduce the impact of dam failure.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare no conflicts of interest regarding the publication of this paper.
Acknowledgments
The financial supports from the National Key Research and Development Plan of China (grant nos. 2018YFC0604605 and 2018YFE0123000) and the Science and Technology Innovation Fund of BGRIMM Technology Group (grant no. JTKJ1812) are gratefully acknowledged.