Geofluids
Volume 2019 (2019), Article ID 6127174, 11 pages
https://doi.org/10.1155/2019/6127174

Research Article

Numerical Simulation of Ore Particle Flow Behaviour through a Single Drawpoint under the Influence of a Flexible Barrier

Qingfa Chen,1,2 Shikang Qin,¹ and Qinglin Chen¹

¹College of Resources, Environment and Materials, Guangxi University, Nanning, Guangxi 530004, China
²State Key Laboratory of Geomechanics and Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences, Wuhan 430071, China

Correspondence should be addressed to Qingfa Chen; chenqf@gxu.edu.cn

Received 11 September 2018; Revised 6 November 2018; Accepted 21 November 2018; Published 5 February 2019

Academic Editor: Mohammad Sarmadivaleh

Copyright © 2019 Qingfa Chen et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

A model of the mass draw and synchronous filling no-top-pillar shrinkage stoping method was created based on the particle element theory and PFC2D software. A numerical simulation of the ore particle flow behaviour through a single drawpoint under the influence of a flexible barrier was performed; following which, a visualization of the change process in the marked layers, draw column, flexible barrier, and cavity was obtained. Finally, the results demonstrated the following. (1) Before the ore particles at the highest horizon were drawn, the morphology of the draw column had an ellipsoidal (whole, closed, and approximate) shape, but after the ore particles were drawn, the morphology of the draw column had a gyroscopic shape (i.e., the shape of a spinning top). (2) The height of the draw column increased rapidly and exponentially in the early period of the draw, presented a linear growth tendency during the medium period with an increase in the accumulated mass drawn, and remained at a certain value in the later period. (3) The morphology of the flexible barrier presented a Gaussian distribution on the whole, but its bottom had an analogous parabola morphology, particularly in the late stage. (4) The cavity began to form with the ore particles at the highest horizon drawn; the cavity development presented an evolutionary process from the microscopic to macroscopic regime, and the morphology of the cavity was crescent in the early stage and triangular in the late stage.

1. Introduction

Under the joint guidance of green, no-waste, and synergetic mining [1–4], the authors proposed the idea of “synchronous filling” in mining technology in 2010 and invented the mass draw and synchronous filling with no-top-pillar shrinkage stoping method (i.e., synchronous filling mining), as shown in Figure 1 [5]. This method is based on the traditional shrinkage method; before drawing the ore, workers can stand on the ore pile (6 in Figure 1) to lay the flexible barrier. In the process of drawing, the filling materials are filled into the stope room through the return airway (1 in Figure 1) in a timely fashion, and by means of vibration and gravity, the filling materials and ore particles sink synchronously and uniformly [6]. The ore particles are drawn under the influence of a flexible barrier, which does not come into direct contact with the waste rocks. Accordingly, there is no ore loss or dilution in the ore particles drawn. The flow behaviour of the ore particles in the new mining method is a breakthrough over traditional draw theories. Hence, the study of the flow behaviour of ore particles through a single drawpoint under the influence of a flexible barrier is meaningful. Currently, there are three major types of classical draw theories: the stochastic medium draw theory, the ellipsoid draw theory, and the similar ellipsoid draw theory [7]. The stochastic medium draw theory [8] simplifies granular media to a continuous-flow stochastic medium and studies the moving process of granular media using a probability method. In the ellipsoid draw theory [9], the mathematical model of the draw column is an ellipsoid equation, and the other regularity phenomena are expressed according to the basic properties of the draw column. The similar ellipsoid draw theory [10] held that the draw column is a similar ellipsoid and established the ideal function of the draw theory of a similar ellipsoid, which improved the ellipsoid draw theory.

Figure 1: Schematic diagram of the mass draw and synchronous filling no-top-pillar shrinkage stoping method: (1) return airway, (2) top pillar, (3) raise, (4) cross-heading, (5) barrier pillar, (6) remaining ore, (7) drift pillar, (8) drawpoint, (9) haulage drive, (10) unmined ore, (11) filling material, and (12) barrier.

Numerical simulation technology has been widely used in draw studies [11–13]. The draw test using PFC software is based on discrete element theory [14] and can essentially analyse and describe the flow law of the caved ore from a microscopic perspective. It intuitively shows the moving, recovery, and residue of ore particles, as well as the rock mixing process. Mining researchers around the world have made numerous achievements in the study of drawing using PFC software. For example, Pierce and Bridgwater [15] studied the seepage problem of fine particles in the draw process using the PFC3D program and proposed a seepage velocity equation based on the experimental results. Based on the numerical simulation results of this software, Lorig and Cundall [16] developed a rapid emulator based on PFC3D REBOP that could quickly simulate the draw process. Liu et al. [17] performed an optimization study on the structure parameters of no-bottom-pillar stoping in the Daye iron mine by using PFC2D. Sun et al. [18] studied the flow characteristics of caved ore rock under complex boundary conditions by constructing a model using PFC3D.

To this end, this paper studied the flow behaviour of ore particles in a single drawpoint under the influence of a flexible barrier by constructing a model of the synchronous filling mining method based on the PFC2D program for a microscopic discrete element method.

2. Model Construction and Parameter Selection

According to the dimensions of the structure outline of the physical experimental model in the synchronous filling mining method, the dimensions of the draw model in this 2D numerical simulation were , and the spaces between drawpoints were 24 cm. The wall structure of the model is shown in Figure 2.

Figure 2: Wall structures of the model.

From Figure 2, the entire structure is composed of 23 walls. Wall numbers 3, 6, 9, 12, 15, 18, and 21 represent the drawpoints in the bottom. Wall numbers 1 and 23 represent the left and right wall of the model, respectively.

To simulate the flow behaviour of ore particles and the evolution law of the barrier associated with a single drawpoint in the synchronous filling mining method, the rolling resistance linear model was used for the contacting of particles, the gravity accumulation method was used for the generation of particles, and the cohering of particles is mainly induced by gravity and boundary in the model. The entire draw process could be divided into the following six stages:

(1)In the model, particles were randomly generated between 8 and 130 cm; a gravitational acceleration of 9.81 m·s⁻² was applied to the particles. The walls and particles were endowed with mesomechanical parameters, as shown in Table 1. These parameters were used to adjust the model to make the particles compacted(2)Particles above 128 cm were deleted, and the mesomechanical parameters of the particles were reset according to the actual situation. Because the actual ore particles are angular, the rolling resistance coefficient was added, as shown in Table 2(3)A layer of particles (0.015 m radius) was formed on the surface of the particles, and the radius of these particles was determined by a thickness of 0.03 m of the barrier used in the physical tests [19, 20]. The parallel bond model was used to simulate the flexible barrier, which was endowed with the mesomechanical parameters shown in Table 3. The value of the elastic modulus of the parallel bond of the flexible barrier is , which is determined according to the trial and error method, and the goal of trial and error is to find a suitable elastic modulus value that can prevent the rupture of the flexible barrier in the simulation process. By inputting the parameters in Table 3 into the software, the particles can be bonded together automatically(4)The particles were grouped at regular intervals in the vertical direction by PFC [21] commands to realize the visualization of the marked particle layer, as shown in Figure 3(5)After the bottom wall (no. 12) was deleted, the particles were drawn continuously, and then the draw process began. The loop statements were compiled with the use of the FISH language in PFC2D [21] to realize the entire draw process. The function of the loop statement is as follows: (1) to simulate the synchronous filling draw, the filling waste particles were properly generated on the surface of the particles after a certain number of particles were drawn, and the mesomechanical parameters of the waste particles are shown in Table 4, and (2) the ID numbers of the particles drawn were recorded until the particles were not drawn in the model(6)The initial balance position of each drawn particle could be obtained through the coordinate value (x, y) of the particle and the ID number of the particle drawn, and the region formed by these particles was the draw column. Through the above processes, the visualization of the morphology of the draw column could be achieved, and the concrete position of the draw column in the stope could be described

Table 1: Mesomechanical parameters of the walls and initial particles.

Wall
Initial particles

Normal stiffness (N/m)
Shear stiffness (N/m)
Friction coefficient
Normal stiffness (N/m)
Shear stiffness (N/m)
Friction coefficient
Density (kg/m³)
Radius (m)

0.5

0.3
2800
0.008

Table 2: Mesomechanical parameters of the particles.

Normal stiffness (N/m)
Shear stiffness (N/m)
Friction coefficient
Density (kg/m³)
Rolling resistance coefficient
Radius (m)

0.5
2800
0.5
0.008

Table 3: Mesomechanical parameters of the barrier.

Normal stiffness (N/m)
Shear stiffness (N/m)
Normal stiffness of the parallel bond (N/m)
Shear stiffness of the parallel bond (N/m)
Density (kg/m³)
Friction coefficient
Elastic modulus of the parallel bond (N/m²)
Radius (m)

2000
0.4

0.0015

Figure 3: Layers of marked particles in the model.

Table 4: Mesomechanical parameters of the waste particles.

Normal stiffness (N/m)
Shear stiffness (N/m)
Friction coefficient
Density (kg/m³)
Rolling resistance coefficient
Radius (m)

0.5
2800
0.5
0.008

3. Analysis of Numerical Simulation Result

3.1. Test Phenomena

In the draw simulation test associated with a single drawpoint, the particles were continuously drawn, the barrier gradually sank, and the interface between the particles and waste rocks was kept flat by synchronously filling waste rocks. After a certain amount of ore was drawn, the relevant data, such as the number and location of the marked particles drawn, were recorded. The flow pattern of particles is shown in Figure 4.

(a) Fifth draw

(b) Tenth draw

(c) Fifteenth draw

(d) Twentieth draw

(e) Twenty-fifth draw

(f) Last draw

Figure 4: Flow pattern of the particles.

After opening the drawpoint, the particles near the drawpoint began to be drawn, and the range expanded with an increasing number of particles drawn. The marked particles of each layer sank successively from the bottom to top, and the morphology of which presented a Gaussian distribution on the whole. The morphology of the barrier after sinking was similar to that of the marked particles, which also presented a Gaussian distribution pattern. However, the morphology of the barrier at the bottom was relatively smooth, and a cavity formed near the bottom. The morphology of the cavity became more obvious with an increase in the particles drawn. The distance between the marked particles in each layer on the generatrix of the draw cone increased gradually, the interval between the marked layers decreased gradually, and the particles undrawn at the top of the drawpoint were disordered.

According to the traditional draw theory, the marked particles at each layer show a cone-shaped drop in the draw process with a single drawpoint. The marked particles in each layer begin to form a breaking cone when the marked particles in the generatrix reach the drawpoint. The breaking funnel in each layer is not formed until that of the marked particles at the highest horizon has formed. If the draw continues, waste rock will mix in. However, if a barrier is placed on the top of the ore, the waste rock can be isolated and the draw can continue. To explore the flow behaviour of particles in a drawpoint under the influence of a flexible barrier, the morphology of the draw column, the relation between the height of the draw column and the accumulated mass drawn, the moving law of the flexible barrier, and the evolution law of the cavity were analysed.

3.2. Morphology Evolution Law of the Draw Column

Based on the ID number of each particle drawn in the numerical test and the coordinate values of the particles in the initial equilibrium state, the initial equilibrium position of each particle drawn can be obtained by grouping and restoring the particles in the initial experimental model, and the zone formed by these particles is termed a draw column or isolated extraction zone. The morphology evolution process of the draw column in a single drawpoint under the influence of the flexible barrier is shown in Figure 5. The boundary of the draw column was connected smoothly in several typical morphologies, and the diagrammatic plan of the morphology evolution process of the draw column was obtained, as shown in Figure 6.

Figure 5: Morphology evolution process of the draw column in the test.

Figure 6: Diagrammatic plan of the morphology evolution process of the draw column in the test: (a) transitional ellipsoid, (b) limit ellipsoid, (c and d) transitional gyrostat, and (e) limit gyrostat.

The morphology of the draw column was regarded as having an approximately ellipsoidal shape before the particles in the highest horizon were drawn, and the draw column increased gradually to the highest level, as shown in Figures 6(a) and 6(b). The basic rule has not been changed by the existence of the flexible barrier, which is still in accordance with the traditional draw theory. After the particles at the highest level were drawn, the upper part of the draw column became a new curve because of the sliding effect of the barrier and the morphology was no longer a part of the ellipsoidal shape. The middle part of the draw column had an approximately ellipsoidal shape because a cavity existed at the bottom of the barrier, the lower part of the draw column had an extended approximately ellipsoidal shape because it was less affected by the barrier, and the entire draw column had a gyroscopic shape, as shown in Figures 6(c), 6(d), and 6(e).

(1)According to the draw column equation of a nonpoint source draw in the random medium draw theory [22, 23], the morphology of the draw column before the particles in the top level were drawn can be determined:

where is the transverse length of the marked particles in the draw column, is the horizontal height of the marked particles in the draw column, is the height of the draw column, is the height influenced by the drawpoint, and and are the flow parameters of granular media (ore particles) in the stochastic medium draw theory [22, 23]. , where is the diameter of the drawpoint and is the internal friction angle    The plane coordinate system was established where the origin is the lowest point of the draw column. The y-axis goes through the origin and moves straight up, and the x-axis goes through the origin and moves horizontally to the right. When the height of the draw column is 123 cm in the extraction conditions influenced by the flexible barrier, the test data in the boundary of the draw column is shown in Table 5.    By substituting the data in Table 5 into equation (1), the following results can be obtained: , , and the correlation coefficient is 0.996.    Therefore, the equation of the draw column before the particles at the top level were drawn is

(2)After the particles in the highest level were drawn, the morphology of the draw column changed markedly and generally presented a gyroscopic shape, which was quite different from the three classical draw theories, in which the shape of the draw column is ellipsoidal or approximately ellipsoidal. The particles in the upper part were drawn in advance due to the friction action of the barrier, the particles in the middle can easily fall from the boundary of the cavity to the bottom of the cavity and be drawn in advance because of the existence of the cavity, and the flow behaviour of the particles in the lower part remained unchanged. The morphology of the draw column in the upper part was similar to the morphology of the corresponding barrier curve, the morphology in the middle was similar to that in the upper part which is ellipsoidal in shape, and the morphology in the lower part was the extension of the ellipsoidal shape, but it was a short straight line (the inclination of the line is the natural repose angle) in the last draw due to the existence of the cavity. Therefore, the overall gyroscopic shape was composed of an exponential curve in the upper part, an oblique ellipsoidal shape in the middle, and an extension of the ellipsoidal shape in the lower part

Table 5: The test data in the boundary of the draw column.

0
10.68
16.42
16.91
18.07

126.54
114.24
92.18
81.94
56.27

The gyroscopic shape is symmetric with respect to the y-axis, and the final shape of the draw column can be described by the data on the right side. After fitting the three segments of the shape of the draw column, its approximate mathematical expression was obtained:

where is the natural repose angle and . The other symbols have the same meaning as previously described.

The comparison between the curve of the draw column on the right side and the fitting curve is shown in Figure 7.

Figure 7: Comparison graph of the curve of the draw column on the right side and fitting curve at the end of draw.

3.3. Relationship between the Height of the Draw Column and the Accumulated Mass Drawn

In the first 15 draw processes, the particles at the highest level were not drawn. Therefore, the shape of the draw column was ellipsoidal, and the height of the draw column increased with the increase in the accumulated mass drawn, which is in accordance with the traditional draw theory. After the number of the draws reached 15, the particles at the highest level were drawn and the shape of the draw column was transformed into the gyroscopic shape. Meanwhile, the height of the draw column no longer increased with the increase in the accumulated mass drawn (it became a certain value), which is quite different from the traditional draw theory. However, the draw did not stop, and many particles were drawn from the drawpoint.

The two indexes, the height of the draw column and the accumulated mass drawn, were statistically analysed, and the statistical results are shown in Table 6.

Table 6: Statistical results of the height of the draw column and the accumulated mass drawn.

Draw number
Height of the draw column (cm)
Accumulated mass drawn (kg)
Draw number
Height of the draw column (cm)
Accumulated mass drawn (kg)

1
23
84
15
127
1299

2
40
195
16
128
1403

3
53
268
17
128
1499

4
60
308
18
128
1556

5
68
398
19
128
1603

6
77
477
20
128
1659

7
82
554
21
128
1761

8
87
654
22
128
1873

9
95
749
23
128
1937

10
101
853
24
128
2039

11
108
930
25
128
2118

12
113
1024
26
128
2204

13
117
1107
27
128
2226

14
124
1212
28
128
2228

3.3.1. Analysis of the Model Reliability

Castro et al. [24] performed a 3D physical draw experiment presenting the flow characteristics of ore and rock in the caving mining method by a medium of gravel, which has been the largest scale to date. The relationship between the height of the draw column () and the accumulated mass drawn () from the experiment satisfied the following equation:

where and represent, respectively, the height of the draw column and accumulated mass drawn at which the height increases exponentially with the mass and is the final growth velocity when the height grows linearly with the accumulated mass drawn.

To analyse the suitability and reliability of the model, based on the Levenberg-Marquard (L-M) algorithm [25], the relationship between the height of the draw column and the accumulated mass drawn in the first 15 draw processes was fitted nonlinearly in this simulation test. Comparisons between the theoretical curves and the experimental data are shown in Figure 8. The fitting results of the coefficient are shown in Table 7. The (goodness of fit coefficient) in Table 7 is close to 1, indicating that equation (4) complied with the statistical data in the test. Hence, this numerical simulation can truly reflect the draw rules.

Figure 8: Comparison between the theoretical curve and experimental data.

Table 7: Fitting results of the coefficient in equation (4).

56.257
233.588

0.998

3.3.2. Relationship between the Height of the Draw Column and the Accumulated Mass Drawn in the Numerical Test

According to Figure 8 and the fitting results, in the initial stage of draw, the height of the draw column increased exponentially with the accumulated mass drawn; subsequently, with an increase in the accumulated mass drawn, the growth velocity gradually decreased. After increasing to a certain extent, the height of the draw column increased linearly with an increase in the accumulated mass drawn. However, Figure 8 only shows the data from the first 15 draws in line with the traditional ore draw theory and cannot fully explain the relationship between the height of the draw column and the accumulated mass drawn in the entire draw process. Therefore, the data in Table 6 need to be further analysed. The graph of the data in Table 6 is shown in Figure 9.

Figure 9: Relationship between the height of the draw column and the accumulated mass drawn in the numerical test.

Figure 9 describes the relationship between the height of the draw column and the accumulated mass drawn in the numerical test of a single drawpoint. In the figure, the curve before the demarcation line () is the range of the ellipsoid, and the straight line between the demarcation line and the finish line () is the range of the gyrostat. Further combining with Figure 8, under the influence of the flexible barrier, the overall variation law of the height of the draw column with the accumulated mass drawn in a single drawpoint is as follows: before the accumulated mass drawn reached 398 kg, the height of the draw column increased exponentially; when the accumulated mass drawn was between 398 and 1403 kg, the height of draw column increased linearly with the accumulated mass drawn; and after the accumulated mass drawn reached 1403 kg, the height of the draw column did not increase with the increase in the accumulated mass drawn but presented a horizontal linear relation.

3.4. Morphology Evolution Law of the Interface of the Barrier

In the process of performing a numerical simulation on the technology of the mass drawn and synchronous filling, the section data of the barrier interface were recorded after a certain mass of particles had been drawn. The dynamic subsidence graphs of the barrier interface are shown in Figure 10.

(a) Fifth draw

(b) Tenth draw

(c) Fifteenth draw

(d) Twentieth draw

(e) Twenty-fifth draw

(f) Last draw

Figure 10: Dynamic subsidence graphs of the barrier interface.

In Figure 10, the barrier was horizontal before the draw started. After opening the drawpoint, the barrier decreased and bent gradually as the particles in the model were constantly drawn under the combined action of the load caused by backfilling waste particles and the flowing field of particles. When the accumulated mass drawn reached a certain value, an obvious cavity appeared in the lower part of the barrier, and the volume of the cavity gradually expanded. On the profile of the model, the shape of the barrier interface presented a Gaussian distribution curve on the whole. In the draw process, the bottom of the barrier and the particles were separated due to the existence of the cavity, and the bottom of the barrier presented an analogous parabolic morphology under the interaction of the overlying particle pressure and tension from the barrier on both sides, particularly at the later stage.

3.5. Formation Mechanism and Morphology Evolution Law of the Cavity

3.5.1. Formation Mechanism of the Cavity

The formation of the cavity is due to the sinking rate of the barrier being lower than the rate of the particles below the barrier, and the load of the filling waste rocks is not enough to make the barrier sink and make full contact with the particles below the barrier. Therefore, it is easy to separate the barrier from the particles. The sinking rate of the barrier during the draw process is less than that of the particles below the barrier. Therefore, the initial position of the cavity can be explained by the velocity field, and the following assumptions are made:

(1)The intersection of the generatrix of the depression cone and the initial horizontal line of the barrier is supposed to be a particle(2)The initial movement stage of the barrier does not have any influence on the flow behaviour of particles; thus, the particles are in direct contact with the waste rocks in the initial model(3)The particles are nonexpansive, and the density of each point remains unchanged(4)The impulsive load caused by the backfill of the filling material is not considered.

According to the velocity equation on the flow axis of the unexpanded material in random medium draw theory [26], the velocity on the flow axis in the model can be expressed as

where is the vertical falling speed of the particle, is the accumulated mass drawn per unit time, and is the height of the draw column.

Suppose the height of the subsurface of the barrier is and the spacing between the barrier and the particles immediately below the subsurface of the barrier is . The height of the particles adjacent to the subsurface of the barrier is therefore . Based on the assumption conditions and equation (5), the rate of the barrier and the rate of the particles adjacent to the subsurface of barrier can be obtained:

Equation (6) is divided by equation (7):

Therefore, when the particles at the top level first begin to move, the rate of the hypothetical particle at the intersection of the generatrix of the depression cone and the initial horizontal line of the barrier is less than that of the particles adjacent to the subsurface of the barrier. Accordingly, a small void begins to appear between the barrier and the particles. With the gradual subsidence of the barrier, its restrictive effect gradually became obvious, and the speed difference between the two gradually increases, which is the reason for the macroscopic cavity in the later stage of the experiment.

3.5.2. Morphology Evolution Law of the Cavity

When the loose range of particles reaches the top of the particles in the model, the barrier begins to sink. Since the sinking rate of the barrier is less than the rate of the particles adjacent to the subsurface of the barrier, a cavity appears at the bottom of the barrier. The development of the cavity is an evolution process from the microscopic to macroscopic regime, where the shape of the cavity is crescent in the early stage and triangular in the late stage. When the barrier drops to a certain depth, a clear cavity will appear at the bottom of the barrier, as shown in Figure 11.

(a)

(b)

Figure 11: Morphology evolution of the cavity.

Due to the existence of the barrier, particles on both sides of the cavity will continuously leak out in an arch morphology, while waste rock particles will not be drawn from the model. Therefore, the existence of the flexible barrier can reduce the loss and dilution of the ore and allow the pure ore under the barrier to be fully drawn.

4. Conclusions

(1)Numerical simulation test studies on the mass draw and synchronous filling shrinkage method were performed. The flow behaviour of ore particles in this method was recognized. The evolution laws of the draw column, the barrier, and the cavity were obtained(2)The morphology of the draw column has an ellipsoidal (whole, closed, and approximate) shape before the ore particles at the highest horizon have been drawn, but it transforms to have a gyroscopic morphology after the ore particles above are drawn(3)Before the ore particles are drawn at the top level, the height of the draw column increases with an increase in the accumulated mass drawn (it grows exponentially before it grows linearly). After the ore particles are drawn at the top level, the height of the draw column does not increase with the increase in the accumulated mass drawn (it is a certain value)(4)The barrier has a general Gaussian morphology, but the bottom of the layer presents a morphology of an analogous parabola(5)The formation of the cavity is due to the sinking rate of the barrier being lower than the rate of the ore particles below the barrier, and the load of the filling waste rocks is not enough to make the barrier sink and make full contact with the ore particles below the barrier. The development of the cavity presented an evolutionary process from the microscopic to macroscopic regime

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 that they have no conflicts of interest.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant No. 51464005) and the Open Research Fund of State Key Laboratory of Geomechanics and Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences (Grant No. Z016015).

References

G. Liying, R. Zhiyuan, and L. Yansui, “The causes of land landscape changes in semi-arid area of Northwest China: a case study of Yulin City,” Journal of Geographical Sciences, vol. 16, no. 2, pp. 192–198, 2006.

L. Ru-yin and Z. Xiao-ting, “Negative entropy mechanism of the circular economy development countermeasures in mining area,” Procedia Earth and Planetary Science, vol. 1, no. 1, pp. 1678–1685, 2009.

Z. Bian, H. I. Inyang, J. L. Daniels, F. Otto, and S. Struthers, “Environmental issues from coal mining and their solutions,” International Journal of Mining Science and Technology, vol. 20, no. 2, pp. 215–223, 2010.

B. Milanez and J. A. Puppim de Oliveira, “Innovation for sustainable development in artisanal mining: advances in a cluster of opal mining in Brazil,” Resources Policy, vol. 38, no. 4, pp. 427–434, 2013.

Q. F. Chen and Z. X. Wu, A Large Number of Ore Drawing Synchronous Filling No-Top-Pillar Shrinkage Stopping Method, CN, 2010.

Q. F. Chen and Q. L. Chen, “Synchronous filling mining technology idea and a kind of representative mining method,” China Mining Magazine, vol. 24, pp. 86–88, 2015.

G. M. Malakhov, “Theoretical principles of ore drawing and factors influencing index of yield,” Powder Technology, vol. 3, no. 1, pp. 364–366, 1969-1970.

J. Litwiniszyn, “Application of the equation of stochastic processes to mechanics of loose bodies,” Archives of Mechanics, vol. 8, no. 4, pp. 393–411, 1956.

T. M. Malakhov, Drawing of Caving Ore Block, Metallurgical Industry Press, Beijing, 1958.

R. F. Li, “Actual equation of the theory of ore drawing of similar ellipsoid,” Nonferrous Metals(Mining), vol. 45, pp. 36–42, 1994.

R. L. Castro, F. Gonzalez, and E. Arancibia, “Development of a gravity flow numerical model for the evaluation of drawpoint spacing for block/panel caving,” Journal of the Southern African Institute of Mining and Metallurgy, vol. 109, pp. 393–400, 2009.

G. P. Chitombo, “Cave mining: 16 years after Laubscher’s 1994 paper ‘Cave mining - state of the art’,” Mining Technology, vol. 119, no. 3, pp. 132–141, 2010.

H. Sun, A. B. Jin, Y. T. Gao, and X. Q. Meng, “Flow characteristics of caved ore and rock in the multiple draw-point condition,” Chinese Journal of Engineering, vol. 37, p. 1259, 2015.

P. A. Cundall and O. D. L. Strack, “A discrete numerical model for granular assemblies,” Géotechnique, vol. 29, no. 1, pp. 47–65, 1979.

M. Pierce and J. Bridgwater, PFC3D Modelling of Inter-Particle Percolation in Caved Rock under Draw. Numerical Modeling in Micro Mechanics via Particle Methods, Taylor & Francis Group, London, 2004.

L. J. Lorig and P. A. Cundall, “A rapid gravity flow simulator,” Tech. Rep., JKMRC and Itasca Consulting Group Inc, Brisbane, 2000.

Z. N. Liu, L. F. Mei, and W. D. Song, “Study on stope structure parameter optimization of sublevel caving without sill pillar based on numeral simulation with PFC,” Mining Research & Development, vol. 28, p. 3, 2008.

H. Sun, A. B. Jin, Y. Q. Gao, Y. Zhou, and Z. Yang, “Flow characteristics of caved ore and rock under complex boundary conditions,” Journal of Central South University (Science and Technology), vol. 37, pp. 1251–1259, 2015.

Q. F. Chen, Q. L. Chen, J. Y. Zhong, D. P. Chen, S. H. Li, and W. J. Niu, “Evolution law of interface morphology of flexible isolation layer under ore drawing from single funnel,” The Chinese Journal of Nonferrous Metals, vol. 26, pp. 1332–1338, 2016.

Q. F. Chen, Q. L. Chen, J. Y. Zhong, Y. D. Wang, S. H. Li, and W. J. Niu, “Flow pattern of granular ore rock in a single funnel under a flexible isolation layer,” Chinese Journal of Engineering, vol. 38, 898 pages, 2016.

IC Group, PFC3D (Particle Flow Code in 3 dimensions) FISH in PFC3D, Itasca Consulting Group Inc, Minnesota, 2008.

G. Q. Tao, F. Y. Ren, and Z. D. Liu, “Research on improvement of the stochastic medium theory for ore drawing,” Journal of Mining and Safety Engineering, vol. 27, pp. 239–243, 2010.

G. Q. Tao, S. J. Yang, and F. Y. Ren, “Test on granules flow parameters of ore drawing with stochastic medium theory,” Chinese Journal of Rock Mechanics and Engineering, vol. 28, pp. 3464–3470, 2009.

R. Castro, R. Trueman, and A. Halim, “A study of isolated draw zones in block caving mines by means of a large 3D physical model,” International Journal of Rock Mechanics and Mining Sciences, vol. 44, no. 6, pp. 860–870, 2007.

L. Chang, J. Xie, S. Wang, and J. Xiao, “The simulation and prediction of groundwater regime based on the optimized algorithm of L-M neural network,” Groundwater, vol. 27, pp. 380–383, 2005.

F. Y. Ren, Study and Application of Stochastic Medium Theory in Ore-Drawing, Metallurgical Industry Press, Beijing, 1994.