*Article* **Numerical Investigation on the Impact of Tailings Slurry on Catch Dams Built at the Downstream of a Breached Tailings Pond**

**Shitong Zhou 1,2 and Li Li 2,\***


**Abstract:** Tailings storage facilities (TSFs) are known as a time-bomb. The numerous failures of TSFs and the heavy catastrophic consequences associated with each failure of TSFs indicate that preventing measures are necessary for existing TSFs. One of the preventing measures is to construct catch dams along the downstream near TSFs. The design of catch dams requires a good understanding of the dynamic interaction between the tailings slurry flow and the catch dams. There are, however, very few studies on this aspect. In this study, a numerical code, named LS-DYNA, that is based on a combination of smoothed particle hydrodynamics and a finite element method, was used. The numerical modeling shows that the tailings slurry flow can generally be divided into four stages. In terms of stability analysis, a catch dam should be built either very close to or very far from the TSF. When the catch dam with an upstream slope of a very small inclination angle is too close to the tailings pond, it can be necessary to build a very high catch dam or a secondary catch dam. As the impacting force can increase and decrease with the fluctuations back-and-forth of the tailings slurry flow, the ideal inclination angle of the upstream slope of the catch dam is between 30◦ and 37.5◦, while the construction of a catch dam with a vertical upstream slope should be avoided. However, a catch dam with steeper upstream slopes seems to be more efficient in intercepting tailings flow and allowing the people downstream to have more time for evacuation. All these aspects need to be considered to optimize the design of catch dams.

**Keywords:** tailings dam; impacting force; kinetic energy; numerical simulation; smoothed particle hydrodynamics (SPH)

#### **1. Introduction**

Mines produce large volumes of mine waste each year in terms of tailings and waste rocks [1]. A small portion of these materials can be sent back to underground voids as backfill [2–18]. The recycling and reuse of these materials for in- and out-mine sites are also increasingly seen [19]. All these practices can reduce the amount of mine waste to be disposed of on the ground surface. However, due to the swelling of rocks by blasting and grinding, the large amount of tailings has to be sent by pipes to a place on the ground surface contained by dams, called a tailings storage facility (TSF), which then becomes a new perpetuity infrastructure in the environment [20–25].

Tailings are a byproduct of mill processing, containing a large portion of fine particles smaller than 80 μm and water [26]. Upon disposition in a TSF, tailings remain saturated for a long time. This is particularly true in the regions where the water balance is positive. When a TSF fails, a large volume of tailings can be instantaneously released. Their flow usually transforms into a mud flood, causing significant damage to the environment and infrastructures, and even loss of lives [27]. In 1965 in Chile, the town of El Cobre

**Citation:** Zhou, S.; Li, L. Numerical Investigation on the Impact of Tailings Slurry on Catch Dams Built at the Downstream of a Breached Tailings Pond. *Processes* **2022**, *10*, 898. https://doi.org/10.3390/pr10050898

Academic Editor: Farhad Ein-Mozaffari

Received: 13 February 2022 Accepted: 29 April 2022 Published: 2 May 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

was destroyed and more than 200 people were killed due to the failure of an upstream TSF [28]. Since then, a lot of work has been done on the analysis of failure mechanisms of TSFs [29–33]. Despite considerable progress in the studies, failures of TSF continue to take place worldwide. In February 1994, in South Africa, the village of Merriespruit was destroyed and 17 people were killed due to the failure of an upstream TSF. On that occasion, the TSF was, for the first time, qualified by a judge as a time-bomb [34]. In China, the failure of a tailings dam at Shiqiaozi, Haicheng, Liaoning Province, on 25 November 2007 resulted in 10 people being killed, 3 reported missing, and 17 injured. Another case of a TSF failure at Taoshi, Linfen City, Xiangfen County, Shanxi Province, on 8 September 2008 involved 277 people being killed and 33 injured [28]. In 2014, in Canada, the failure of a TSF of the Mount Polley mine resulted in the release of 23.6 million cubic meters of tailings into the environment [35,36]. In November 2015, at the Samarco mine in Brazil, the failure of a TSF resulted in the death of 17 people and a flow of tailings sludge of about 800 km up to the Atlantic coast [28,37,38]. In January 2019, near Brumadinho in Brazil, 248 people were killed and 22 were reported missing due to the failure of a TSF [28,39]. These endless failures of TSFs worldwide indicate that big efforts are still needed to prevent more casualties issued from TSF failures.

The instability of a TSF can be related to a number of influencing factors, varying from technical and engineering to maintenance and management [27,29,30,34,38–41]. The technical and engineering aspects affecting the stability of TSFs include the design and construction of TSFs and their foundations. At current stage, it is impossible to entirely eliminate any risk of TSF failure [33,42]. Rather, the probability of TSF failures can even increase with time due to the persistence of existing TSFs, the ceaseless increase of the number of mines and TSFs, and the more and more pronounced climate changes, which usually have a negative impact on the stability of TSFs. Efforts should be continuously made to improve the stability of TSFs [43], install monitoring systems [44], and even forecast the scale of destructions associated with the failure of a TSF [45].

Currently, thickened or filtered tailings are widely considered as the solution of solving the physical instability of TSFs. It is, however, problematic with chemically reactive mine waste because the unsaturated tailings favor the oxidation and generation of acid mine drainage [46–48]. In addition, unsaturated tailings can also be prone to liquefaction, as long as the degree of saturation approaches full saturation [49], and this is probably the case in the regions where the water balance is positive. Considering this as a sure and certain solution to avoid a TSF failure could be irresponsible [21]. Moreover, work is needed to secure existing conventional TSFs. In all cases, risk analysis on the breach of TSFs is necessary [50–65]. In many cases, preventing measures can be necessary to avoid catastrophic consequences in case of a TSF failure. One of the preventing measures is to build catch dams near the downstream of a TSF. Upon failure of the TSF, the catch dams must remain stable to hold at least part of the released tailings slurry, delaying and reducing the flood of tailings. If the catch dam is not properly designed and constructed, it can fail upon the dynamic impact of a tailings slurry flood. This was the case after the breach of a TSF in March 2020 in Yichun, Heilongjiang, China. Ten catch dams had to be built, one after another, along the downstream to stop the tailings flood [66].

The purpose of this study is to obtain a good understanding of the characteristics of the tailings slurry flow after a breach of a TSF, and on the dynamic interaction between the flooding tailings and catch dams. This can ultimately help to obtain a safe and economic design of catch dams to avoid catastrophic consequences after the failure of a TSF.

Over the years, a number of works have been published on the dynamic interaction between water and downstream obstacles [67–72]. One sees also some publications on the dynamic interaction between debris flow and check dams [73–76], or even between granular material flows and obstacles [77]. These works provide a useful insight into the dynamic interaction of water, debris, and granular material with obstacles. However, very few works have been made on the dynamic interaction between tailings slurry and downstream obstacles.

Wang et al. [78] conducted numerical modeling for a specific project to study the flow characteristics of tailings upon the failure of a tailings dam. The effects of downstream dams were briefly presented without any detailed quantitative studies. Zeng et al. [79] evaluated the effect of check dam height against the debris flow of tailings slurry. A similar work has been published by Wang et al. [80], who analyzed the risk associated with the possible failure of tailings dams and the effect of a preventing dam height, built 110 m downstream from the tailings pond. All these works provide some insights on the effects of check dams on the flow of the tailings slurry. However, they were mostly for some specific projects with site-specific conditions. Much more work is needed to obtain a better understanding of the effect of downstream catch dams on the flow of a tailings slurry.

In this study, the dynamic interaction between a tailings slurry flow and a downstream infrastructure is investigated through numerical modeling, using a hybrid numerical approach with the smoothed particle hydrodynamics-finite element method (SPH-FEM). The kinetic energy of the tailings flow and the striking forces on the downstream infrastructure are evaluated. The influences of the position and upstream slope angle of the catch dam on the impact of the catch dam are, for the first time, studied. The results provide some guidelines on the layout and design of catch dams.

#### **2. Numerical Code**

Studying the characteristics of tailings slurry flow needs a numerical code that can simulate the large deformation and displacement with the time of fluid-like materials, while the fluid–structure interaction behavior can be evaluated by the conventional finite-element or finite difference method.

In this study, a numerical code, called LS-DYNA [81], based on a combination of smoothed particle hydrodynamics (SPH) and the finite element method (FEM) was used. SPH is a meshless particle method based on the Lagrangian formulation for fluid dynamics' simulation. It was originally developed by Lucy [82], Gingold, and Monaghan [83] for three-dimensional astrophysical analysis. It is based on the hypothesis that the problem domain can be divided into a set of discrete elements, which possess or can be attributed to a certain material property. With the SPH method, large deformations can be handled without element distortions. It can thus be used to simulate dynamic problems involving large deformations, such as fluid problems with a free surface [84].

The application of the SPH method involves two approximations to represent the simulated space: the kernel and particle [85]. Figure 1 schematically shows a presentation of the SPH kernel function, *W*. On the figure, *r* is the distance between the two particles; *λh* denotes the influence area of the kernel function; *λ* is determined by the kernel function and space dimensions to satisfy the normalization condition.

**Figure 1.** A schematic presentation of the SPH kernel function *W* (reproduced with minor changes from Atif et al. [86], with permission from Copyright Clearance Center, Inc. (CCC) on behalf of Springer-Verlag, Dordrecht, The Netherlands).

In the approximation of the kernel, the spatial distance between the particles is covered by a smoothing length, over which the properties are smoothed by a smoothing kernel function. In particle approximation, the state of a system is represented by a set of particles, which possess individual material properties and move according to the conservation equations. A given particle, *i*, with its given velocity, density, mass, and pressure, interacts with particle *j* within the smoothing length. Considering the conservations of mass, momentum, and energy leads to the following expressions [87]:

$$\frac{d\rho\_i}{dt} = \sum\_{j=1}^{N} m\_j \upsilon\_{ij}^{\beta} \frac{\partial \mathcal{W}\_{ij}}{\partial \mathbf{x}\_i^{\beta}} \tag{1}$$

$$\frac{d\sigma\_i^a}{dt} = \sum\_{j=1}^N m\_j \left( \frac{\sigma\_i^{a\beta}}{\rho\_i^2} + \frac{\sigma\_j^{a\beta}}{\rho\_j^2} + \prod \tau\_{i\bar{\jmath}} \right) \frac{\partial W\_{i\bar{\jmath}}}{\partial \mathbf{x}\_i^{\beta}} \tag{2}$$

$$\frac{de\_i}{dt} = \frac{1}{2} \sum\_{j=1}^{N} m\_j \left( \frac{\sigma\_i^{a\beta}}{\rho\_i^2} + \frac{\sigma\_j^{a\beta}}{\rho\_\flat^2} + \prod\_{i\bar{\eta}} \right) v\_{i\bar{j}}^a \frac{\partial W\_{i\bar{j}}}{\partial x\_i^{\beta}} \tag{3}$$

where *ρ<sup>i</sup>* and *ρ<sup>j</sup>* denote the densities of particles *i* and *j*; *N* denotes the number of particles; *α* and *β* denote the space vectors; *x β <sup>i</sup>* denotes the coordinate component of *x* along *β*; *Wij* denotes the smoothing function between particles *i* and *j*; *v<sup>α</sup> ij* and *v β ij* represent the velocities of particle *i* relative to particle *j* along *α* and *β*; *mj* denotes the mass of particle *j*; *v<sup>α</sup> <sup>i</sup>* and *<sup>v</sup><sup>α</sup> <sup>j</sup>* denote the velocity components of particles *i* and *j* along *α*; *σαβ <sup>i</sup>* and *<sup>σ</sup>αβ <sup>j</sup>* denote the total stress tensors of particles *i* and *j*; ∏*ij* denotes the artificial viscosity; *ei* represents the internal energy of the particle *i*; and *t* represents the time.

While the SPHmethodis advantageousin dealing with fluids andlarge deformation problems, this dynamic meshless technique is much less efficient in dealing with static and small deformation problems of solid bodies than FEM [2,88]. The combination of SPH and FEM methods is thus a natural choice to deal with the dynamic interaction between fluid and infrastructures. This explains the wide application of the SPH and FEM hybrid algorithm in the simulation of the fluid–solid interaction [89].

Figure 2 shows a possible interaction between fluid particles and structure elements. The interaction between fluid and a structure is achieved by a penalty-type 'node-to-surface' contact algorithm. The FEM elements and SPH particles are considered as master and slave parts, respectively. The reaction force between a SPH particle and the contact surface of the structure is then proportional to the depth of the slave part on the surface of the master part. The forces exerted by the SPH particles on the FEM part are then calculated by summing up the forces exerted by each SPH particle on the FEM part.

**Figure 2.** The coupling of SPH and FEM (reproduced with minor changes from Koneshwaran et al. [90], with permission from Copyright Clearance Center, Inc. (CCC) on behalf of Elsevier Science & Technology Journals, Amsterdam, The Netherlands).

In this study, the explicit dynamic solver of LS-DYNA® version R10.0 was employed to perform the numerical modeling based on SPH theories [91]. The applicability (usually called "validation") of the numerical code SPH was tested against the experimental results of a slump test (see details given in Appendix A). The results indicate that the numerical code is capable of simulating the dynamic interaction between a tailings slurry and a structure.

#### **3. Numerical Simulations of Dynamic Interaction between Tailings Flow and A Catch Dam**

#### *3.1. Numerical Models*

Figure 3 shows a numerical model to simulate a largely reduced scale tailings pond. The tailings pond is 0.5 m long, 0.5 m wide, and 0.5 m high. The tailings slurry is the same as that used to test the applicability of the numerical code, presented in Appendix A. It is characterized by a density of 1826.6 kg/m<sup>3</sup> and a dynamic viscosity of 0.43 Pa·s to simulate the tailings slurry having a solids content by mass of 70% of Gao [92]. The catch dam is simulated by a rigid triangular prism 0.2 m high and 0.5 m wide with an upstream slope angle *α* varying from 0◦ to 90◦. The distance between the tailings pond and the catch dam is 0.8 m. The movement of the tailings particles is allowed in all directions, but only inside the open box defined by the back wall of the tailings pond, a horizontal floor and two vertical longitudinal walls. The floor and the three walls are rigid and their movements are not allowed in any direction. The tailings slurry is simulated by particles with SPH, while the rigid catch dam is built with FEM shell elements.

**Figure 3.** The numerical model and boundary conditions of tailings pond break model.

In this study, the flow behavior of the tailings slurry and the impact of the tailings slurry on the catch dam are analyzed by a contact algorithm with the tailings slurry and the catch dam defined as slave and master parts, respectively. The gravitational acceleration is applied along the vertical downward direction. The three vertical walls and the horizontal floor are fixed along their outside surfaces. The catch dam is fixed along its bottom. The simulations were done by considering a duration of 2.5 s, starting from the breach of the tailings pond.

For each numerical model, sensitivity analyses were made to make sure that the distinct elements were dense enough for the tailings and voids simulated by SPH (see Appendix B). These analyses showed that the tailings can be modeled with 357,911 SPH particles and a particle-to-particle spacing of 7 mm.

#### *3.2. Numerical Results*

Before the breach of the tailings pond, the tailings slurry had the highest gravitational potential energy. After the breach of the tailings pond, the tailings slurry flowed. The gravitational potential energy and kinetic energy are interconverted. It is noted that the kinetic energy *K* and gravitational potential energy *P* designate the total kinetic energy and total gravitational potential energy of the whole tailings slurry. At a given time, they can be obtained as follows:

$$K = \sum\_{i=1}^{i=N} \frac{m\_i \upsilon\_i^2}{2} \tag{4}$$

$$P = \sum\_{i=1}^{i=N} \left( m\_i g h\_i \right) \tag{5}$$

where *mi*, *vi*, and *hi* are the mass, velocity, and elevation of SPH particle *i*, respectively; *N* is the total number of SPH particles of the tailings slurry; and *g* is the gravitational acceleration.

Figure 4 shows the flow states of the tailings particles at different times (*T*) starting from the breach of the tailings dam, with the catch dam having an upstream slope angle of *α* = 30◦. At *T* = 0 s, when the tailings dam instantaneously collapses, the potential energy of tailings slurry is the highest and ready to be released. No particle moves and the kinetic energy is zero (Figure 4a). At *T* = 0.24 s, the upper front part of the tailings falls, moves forward and reaches the catch dam with a maximum velocity of 4.40 m/s (Figure 4b). At *T* = 0.5 s, all the upper part of the tailings starts moving. The front line of the tailings flood rushes away from the catch dam at a velocity of 3.31 m/s. Due to this high velocity and the slope inclination of the catch dam's upstream, the movement of the front line is upward, almost parallel to the surface of the upstream slope (Figure 4c). At *T* = 0.76 s, more tailings at the top start to move and the front line of the tailings flood continues rushing away from the catch dam, but falling toward the floor. The maximum velocity of the front line is 3.92 m/s (Figure 4d). At *T* = 1.2 s, all the particles of tailings that can move have moved. It is interesting to note that the tailings slurry near the upstream of the catch dam is not only higher than the catch dam, but also higher than the remnants of the tailings pond. This is probably due to the drag effects by the friction of the upper mobile particles and the suction of front particles. In addition, it should be noted that a back flow of the tailings slurry takes place at the toe of the catch dam's downstream (Figure 4e). This back flow may have a negative impact on the stability of the catch dam. At *T* = 2.5 s, the crest of the tailings flood moves back toward the tailings pond (Figure 4f). This phenomenon is well-known as sloshing [93,94]. Finally, one notes that the tailings particles near the rigid back walls of the tailings pond always show delayed movement, probably due to the friction and viscosity along the tailings–wall interfaces (see Figure 4b,c,f).

Figure 5 shows the variation of the tailings' kinetic energy (Figure 5a) and impacting force (Figure 5b) exerted by the tailings flood on the catch dam's upstream slope surface as a function of time since the breach onset of the tailings dam. One sees that the flow characteristics of the tailings slurry can generally be divided into four stages:

• Stage 1: Initiation. This stage starts at the breach onset of the tailings dam at *T* = 0 s and ends at *T* = 0.38 s when the kinetic energy and impacting forces reach their first peak values. Upon the breach of the tailings dam at *T* = 0 s, the front part of the tailings tends to flow, driven by the potential energy. Once started, the movement of the tailings particles accelerates upon their residual potential energy and the pushing forces exercised by the upper and back particles. Subsequently, the kinetic energy continuously increases with time. At *T* = 0.24 s, the front line of the tailings flood reaches the toe of the catch dam. The upstream slope of the catch dam starts to receive an impacting force. After then, the tailings slurry flows upward over the upstream slope of the catch dam. A part of the kinetic energy has to be transferred into the potential energy required to allow the rising of the tailings particles. This leads to a decreased pace in the increase of kinetic energy until a peak value at *T* = 0.38 s when the impacting force also reaches its first peak value.


**Figure 4.** Flow states of the tailings particles at (**a**) *T* = 0 s, (**b**) *T* = 24 s, (**c**) *T* = 0.5 s, (**d**) *T* = 0.76 s, (**e**) *T* = 1.2 s, and (**f**) *T* = 2.5 s starting from the breach of the tailings dam; simulation with a catch dam having an upstream slope angle of 30◦.

Additional simulations have been done to analyze the effects of the upstream slope inclination angle of the catch dam on the kinetic energy and impacting forces. Figure 6 shows the flow states of the tailings particles when the catch dam's upstream slope starts to be more or less fully covered by the tailings slurry, with the inclination angle *α* varying from 15◦ to 90◦. Firstly, one sees that the tailings flow rises smoothly along the slope when the catch dam's upstream slope has a small inclination angle. Backflow of the tailings front and splashing are observed when the slope inclination becomes abrupt. Secondly, as the distance between the front toe of the tailings and the catch dam's downstream slope toe is kept constant at 0.8 m (see Figure 3), the distance between the front toe of the tailings and the catch dam's upstream slope toe decreases as the inclination angle *α* decreases. The tailings flow thus meets the catch dam and covers the upstream slope earlier with a smaller inclination angle *α*. The catch dam can play its role of resistance earlier, largely reducing the maximum value of velocity. From this point of view, it seems that the catch dam should be constructed as close as possible to the tailings dam. However, when the catch dam with an upstream slope of a very small inclination angle is very close to the tailings pond, a greater quantity of tailings slurry may flow out over the catch dam. Another catch dam along the downstream of the catch dam may be necessary. These results indicate that an optimization process is necessary by considering the position, height, and upstream slope inclination angle of the catch dam, as well as the volume of tailings slurry to be confined.

**Figure 6.** Flow states of the tailings particles when the catch dam's upstream slope starts to be more or less fully covered by tailings slurry, with the inclination angles *α* varying from 15◦ to 90◦: (**a**) 15◦; (**b**) 30◦; (**c**) 45◦; (**d**) 60◦; (**e**) 75◦; (**f**) 90◦.

Figure 7a shows the impacting forces as a function of time for different values of *α*, while Figure 7b presents the peak and trough values of the impacting forces as a function of the upstream slope inclination angle *α* of the catch dam. One sees that the trends of impacting force evolution for different inclination angle values are quite similar (Figure 7a). For a given inclination angle *α*, the impacting force increases to a peak value rapidly after the tailings flow reaches the catch dam. After then, the impacting force decreases to a

minimum value. After a short stage of fluctuations accompanied with slight increases, the impacting force dramatically decreases to a trough value (almost −400 N for the case of *α* = 90◦) and sharply increases to almost the same level before falling. After then, another wave of fluctuations can be observed, but the impacting force generally decreases toward zero. It is interesting to note that these results have a certain similarity with those obtained by the numerical modeling performed by Calvetti et al. [95] on the interaction between the flow of dry granular material and rigid barriers. In the latter case, the impacting force also increases with time until a peak value. After then, the impacting force decreases almost monotonously with time, showing neither a second peak, nor a negative value. This last behavior is a particularity belonging to the flow of tailings slurry associated with the sloshing and possible suction.

**Figure 7.** The Y force exerted on the obstacle along the Y direction: (**a**) time-history of Y force under different inclination angles of the obstacle; (**b**) peak and trough values of Y force as a function of upstream slope inclination angle of catch dam.

By comparing the flow state of tailings and the impacting force evolution, one notes that the impacting forces reach their maximum peak value when the catch dam's upstream slope is almost fully covered by the tailings flow, while the minimum trough value corresponds to the moment of backflow, which tends to pull the upstream slope by suction toward the tailings pond. These results are important for the stability analyses of catch dams. On the one hand, the impacting force should be kept as small as possible in order to avoid pushing the catch dam away, but on the other hand, negative values should be avoided because the upstream slope of the catch dam can be disintegrated by the pulling forces if the catch dam is built with an uncemented or lowly cemented, cohesionless material. The difference between the peak and trough values should also be minimized in order to prevent a cycle loading on the catch dam's upstream slope. From Figure 7b, one sees that the smallest peak value of the impacting force and the smallest difference between the peak and trough values of the impacting force are obtained when the inclination angle *α* is 30◦. The largest peak value of the impacting force and the largest difference are observed when the inclination angle *α* is 90◦. These results tend to indicate that the optimum inclination angle *α* is 30◦, while the construction of the catch dam with a vertical upstream slope represents the highest risk and should be avoided in terms of stability for the catch dam.

Figure 8a shows the evolution of kinetic energy for different upstream slope inclination angles *α*. When the catch dam is absent (i.e., *α* = 0◦), the kinetic energy of the tailings flow increases continuously and reaches a peak value of 256.5 J at *T* = 0.58 s, followed by a gradual decrease due to the frictional resistances along the side walls and floor. These results are very important because they indicate that a catch dam should be built either very close to or very far from the tailings pond. The worst design is to construct the catch dam at a position where the kinetic energy reaches the maximum value and a very strong catch dam is necessary.

**Figure 8.** Kinetic energy of tailings flow under different inclinations of obstacle. (**a**) Time-history of kinetic energy; (**b**) first and second peak of kinetic energy.

From Figure 8a, one sees that two peak values in the kinetic energy can generally be observed when a catch dam is constructed at the downstream of the tailings pond. They can be designated as *K*<sup>1</sup> for the first one and *K*<sup>2</sup> for the second one. The first peak value *K*<sup>1</sup> for the inclination angle *α* = 15◦ is smaller and arrives earlier than the other cases, probably due to the shorter distance between the tailings downstream slope toe and the catch dam's upstream slope toe compared to the other cases. When the inclination angle *α* is larger than 45◦, the magnitude and arrival time of the first peak value *K*<sup>1</sup> become insensitive to further variation in the inclination angle *α*. However, the second peak value *K*<sup>2</sup> arrives later and becomes smaller as the inclination angle *α* increases from 45◦ to 90◦.

Figure 8b shows the variation of the two peak values (*K*<sup>1</sup> and *K*2) as a function of the inclination angle *α*. When *α* = 15◦, both *K*<sup>1</sup> and *K*<sup>2</sup> are small with a relatively large difference between them. At *α* = 90◦, the value of *K*<sup>1</sup> becomes the largest while the value of *K*<sup>2</sup> becomes the smallest (except in the case without a catch dam). This case should be avoided in terms of the stability of the catch dam; this is consistent with the previous analysis of the impacting force. At an inclination angle of *α* = 37.5◦, the magnitudes of *K*<sup>1</sup> and *K*<sup>2</sup> become almost the same even though their values are relatively high. This result is consistent again with the previous analysis of impacting forces as a function of the inclination angle *α*. It seems that the optimal inclination angle *α*, in terms of the stability of the catch dam, ranges from 30◦ to 37.5◦. It is interesting to note that waste rocks typically have a repose angle of 37◦ [96]. This study tends to indicate that waste rocks can be the ideal building materials to construct catch dams. However, detailed analyses should be made in the design of catch dams because waste rocks and other cohesionless materials cannot submit any negative (pulling) forces.

Figure 9 shows the state of the tailings flow at the arrival of the second peak value *K*2. One sees that the front of the tailings flow is closer to the catch dam with a higher value of the inclination angle *α*. Disastrous consequences can be expected if any infrastructures, lives, and equipment are inside the range of the front line. From this figure, one also sees that the delay of the arrival for *K*<sup>2</sup> increases as the catch dam's upstream slope becomes steeper. This is further illustrated in Figure 10, which shows the variation of the arrival time for the second peak value *K*<sup>2</sup> as a function of the inclination angle *α*. When a tailings dam fails, the most important thing is an early warning to the people located along the downstream of the

tailings pond. From this perspective, a larger α seems to be better to intercept the tailings flow and allows the people downstream to have more time for evacuation. However, a steep upstream slope is not the favorite choice for the stability of catch dams. These results show again the complexity in designing a catch dam. Optimization must be done by considering the several influencing factors and possible consequences.

**Figure 9.** Tailings flow forms corresponding to the moment of *K*<sup>2</sup> appearing: (**a**) 15◦; (**b**) 30◦; (**c**) 45◦; (**d**) 60◦; (**e**) 75◦; (**f**) 90◦.

**Figure 10.** Arrival time of *K*<sup>2</sup> as a function of slope inclination angle *α*.

#### **4. Discussion**

In practice, a breach of the tailings dam is usually accompanied with a release of a large quantity of tailings slurry, which transforms into a mud flood and produces extremely destructive consequences on the downstream personnel, infrastructures, and environment. To prevent such undesirable consequences, catch dams can be constructed along the downstream of tailings ponds. The catch dams need to be properly designed. This requires a good understanding on the rheological behavior of the tailings slurry and on the mechanical

responses of the catch dam upon the flooding of tailings slurry. To this end, a series of numerical simulations were performed with LS-DYNA. The numerical results provide some interesting insight on the interaction behavior between the tailings flow and catch dams. The tailings flow can generally be characterized into four stages. Ideal distances between the tailings pond and catch dam should exist. Optimal upstream slope angles of the catch dam have been identified. These numerical results are thus very useful and helpful for the design and construction of catch dams. However, one has to keep in mind that these numerical results were obtained based on some specific case. The study involves several limitations and should be considered as preliminary. Much more work is necessary. For instance, the numerical model used in this study is very small due to the limitation of available computing resources. More numerical simulations with large-dimension models are necessary in the future to analyze the rheological behavior of tailings slurry and the mechanical behavior of catch dams on a large scale.

In this study, the numerical simulations were realized by considering a tailings slurry having a solids content by mass of 70%, a density of 1826.6 kg/m3, and a dynamic viscosity of 0.43 Pa·s [92]. It is well-known that the viscosity and flow characteristics of a tailings slurry depend largely on the solids content, as shown by Chen et al. [65]. More work is necessary to analyze the influence of the solids content on the flow of the tailings slurry and the impact forces on the downstream catch dams.

In practice, the quantity of overflow of the tailings slurry through the catch dam and the maximum distance of the tailings flow, as well as the height of the tailing slurry flood along the two side walls can be of great interest. These are closely related to the height and volume of the tailings pond, the distance between the tailings pond and catch dam, the geometry of the catch dam, and even the downstream topography (side wall inclinations, floor inclination). More work is necessary to better understand the effects of these influencing factors on the tailings flow's characteristics and impact on catch dams in the future.

Another limitation is associated with the simulations of a single catch dam. In practice, the construction of one catch dam may be insufficient to stop the tailings flood, especially when a tailings flood has already started and one has to construct several catch dams somewhere along the downstream far from the tailings ponds [66]. More work is necessary to consider several catch dams positioned at different distances.

Finally, one has to mention that the numerical models have been partly validated against slump tests. Field or laboratory experimental work with tailings ponds and catch dams is necessary to further validate or calibrate the numerical models.

#### **5. Conclusions**

In this work, the flow characteristics of a tailings slurry and its impact force on a catch dam built downstream near the tailings pond after the breach of the tailings dam were analyzed through numerical modeling with a SPH-FEM code, named LS-DYNA. The applicability of the numerical code has been shown through a numerical reproduction of a tailings slurry slump test reported in the literature. Despite the several limitations of the numerical models, several conclusions can be drawn from the numerical modeling performed in this study.

The flow characteristics of the tailings slurry can generally be divided into four stages:


The first analysis shows that the catch dam should be constructed as close as possible to the tailings pond to reduce the maximum velocity of the tailings flow. However, when the catch dam with an upstream slope of a very small inclination angle is too close to the tailings pond, a greater quantity of tailings slurry may flow out over the catch dam. The catch dam needs to be higher or the construction of a secondary catch dam is necessary. Further analysis indicates that a catch dam should be built either very close to or very far from the tailings pond. The worst design is to construct the catch dam at a position where the kinetic energy reaches the maximum value where the construction of a very strong catch dam is necessary.

The numerical results show that the impacting force can increase and decrease with the fluctuations and back-and-forth motion of the tailings slurry flow. A positive and large impacting force may have the effect to push the catch dam away, while a negative force may pull particle material away on the upstream slope of the catch dam. A cycle loading on the catch dam's upstream slope is possible due to the difference between the peak and trough values. From this point of view, the ideal inclination angle *α* is between 30◦ and 37.5◦, while the construction of a catch dam with a vertical upstream slope should be avoided. Waste rocks can be an ideal building material to construct catch dams. However, detailed analyses should be made because waste rocks are cohesionless and cannot submit to any negative (pulling) forces.

While the previous analyses focusing on the stability of catch dams indicate that the abrupt upstream slope of catch dams should be avoided, the arrival time of peak kinetic energy values tends to indicate that a larger value of α seems to be better to intercept the tailings flow and allow the downstream people to have more time for evacuation. Optimization must be done by considering the stability of the catch dam and evacuation time.

**Author Contributions:** S.Z.: conceptualization, numerical modeling, analysis, literature, and writing and editing of the original draft. L.L.: project administration, methodology, supervision, and editing of the original draft. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was financially supported by the China Scholarship Council (202006370203), the Natural Sciences and Engineering Research Council of Canada (RGPIN-2018-06902, ALLRP-566888-21), and industrial partners of the Research Institute on Mines and the Environment (RIME UQAT-Polytechnique; http://rime-irme.ca/). The authors are grateful for their support.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Data are contained within the article.

**Acknowledgments:** The authors gratefully acknowledge the financial support from the NSERC and the industrial partners of the Research Institute on Mining and Environment (RIME). The first author expresses her gratitude for the financial supports from the Central South University (CSU) and the China Scholarship Council (CSC).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Appendix A. Applicability of LS-DYNA to Simulate the Rheological Behavior of Tailings Slurry**

The slump test is widely used to assess the workability or fluidity of slurry or concrete [97,98]. Here, the slump test results of Gao [92] were used to test the applicability of the numerical code LS-DYNA in simulating the rheological behavior of tailings slurry. The slump tests were conducted by following ASTM C-143 with an Abraham cone 300 mm high, 200 mm in the base diameter, and 100 mm in the top diameter. The tailings slurry has a solids content by mass of 70%, a density of 1826.6 kg/m3, and a dynamic viscosity of 0.43 Pa·s.

Figure A1 shows the numerical model to simulate the slump test of Gao [92]. The tailings slurry is treated as an incompressible fluid flow and represented by the material model MAT\_NULL and the equation of state \*EOS\_MURNAGHAN in LS-DYNA [81]. On Figure 3, one sees a meshless infinite rigid wall at the bottom to represent the ground

surface to prevent vertical movement of the fluid passing through it. The friction coefficient between the tailings slurry and the rigid wall is automatically set to be 0.466 [99]. As the numerical model of SPH involves several model controlling parameters [100,101], sensitivity analyses of these model controlling parameters are necessary to ensure stable numerical results. In this study, the sensitivity analyses of model controlling parameters led to a particle spacing of 4 mm (85,980 particles), a smoothing length factor of 1.20, a time step scaling factor of 0.45, and an approximation theory of 15, respectively.

**Figure A1.** The numerical model of the slump test.

Figure A2 shows the variations of the shape and slump height of the tailings slurry with time once the Abraham cone was removed by a lift, obtained by the numerical modeling. As the tailings slurry has a good fluidity, the main slump occurs very quickly within 0.42 s and reaches 273 mm. After then, a further slump becomes extremely slow, with a final slump height of 276 mm. This value is quite close to the test results with a slump height of 286 mm (with a relative error of 3.50%).

Figure A3 shows a comparison between the final shapes of laboratory test and numerical modeling. Both are nearly circular. The base diameter of the slumped slurry can be estimated to be around 350 mm in the photo, whereas the numerical modeling shows a base diameter of 355 mm. The numerical model of LS-DYNA has been proved to be able to simulate the rheological behavior of tailings slurry.

**Figure A2.** The numerical time-history of the slump.

**Figure A3.** Final shape of the tailings slurry after the slump tests: (**a**) a photo of the laboratory test; (**b**) numerical result (Z-coordinate corresponds to the thickness of tailings slurry in meters).

#### **Appendix B. Sensitivity Analyses of SPH Simulations**

In the numerical modeling with SPH, the numerical results can be affected by several controlling model parameters. Particle approximation theory was set as FORM = 15 in this study, according to the suggestion for fluid simulation in LS-DYNA user's manual. Sensitivity analyses of other controlling parameters are necessary to obtain their optimal values and ensure stable numerical results. Figure A4 shows the sensitivity analyses of time step size, the number of particles, and the smoothing length. One sees that the optimal number of SPH particles is 343,000, while the optimal time step scale factor is 0.50. In case of the possible numerical instability caused by a strong splash of the tailings flow, as in the case of declination with 90 degrees, a value of 0.45 was taken for the time step scale factor. Regarding the smoothing length scale factor, the default value of 1.2 was used because the numerical results shown in the figure seemed to not be very sensitive to the variation of this value.

**Figure A4.** Convergence of average kinetic energy of tailings flow under different numerical parameters. (**a**) Mesh independence; (**b**) time step scale factor; (**c**) smoothing length scale factor.

#### **References**

