Abstract

When the upper parts of coal resources are exploited, the mining overburden separated strata may easily connect with the upper unconsolidated aquifers. This may cause the water-sand two-phase fluids in the aquifers to flow underground through the falling broken rock mass. It has been found that it is difficult to eliminate sand after it has entered a mine, which tends to severely restrict the restorations of mines and causes significant economic losses. This paper studied a seepage system dynamic model of water-sand two-phase fluids in broken rock mass based on related two-phase flow theories. In addition, using FORTRAN language, it established a numerical calculation program to study the migration law of water-sand two-phase fluids flowing into broken rock mass and determined the main influencing factors on the seepage characteristic parameters of water-sand two-phase fluids in broken rock mass and the amounts of sand losses.

1. Introduction

Water-sand inrush disasters are common occurrences during the mining of shallow coal seams in China’s western mining areas. Since the coal seam is featured by shallow buried depths and thin base rock, it has been found that, following the exploitation of coal seams, the fallen broken rock mass becomes easily linked with the upper sand aquifers. The water-sand fluids flow underground through the broken rock mass of a mine. This may result in buried equipment on the working surfaces, and even casualties. These negative effects generate major threats to safety, as well as the efficient production of a mine. Figure 1 details the occurring principle of water-sand inrush disaster.

Many researchers have carried out studies which have focused on the water-sand inrush mechanism and have achieved significant results. In summary, the previous research studies were divided into two aspects as follows: first, the development laws of mining fissure or caving zones were examined, for example, the study of the flow channel development laws of water-sand fluids [17]. The second aspect was the study of the water-sand inrush mechanism from the properties of the aquifers, which mainly focused on the examination of the water-rich characteristics of different aquifers and the critical hydraulic slopes of water-sand inrush [8, 9].

In actual situations, with the exception of the water-sand migration channels and the occurrence properties of the aquifers having influences on the water-sand inrush disasters, the migration laws of the water-sand two-phase fluids in broken rock mass are also of great importance for water-sand inrush disasters. These disasters must be studied thoroughly. There have been many previous research studies carried out regarding the seepage in broken rock mass. For example, Martins [10] computed the average seepage velocity of turbulent flow when the Reynolds number was larger than 300 in rubble structures. He studied the impacts of particle size, graduation, and the seepage sections on the seepage in rubble structures and also evaluated stability during seepage. Hansen et al. [11] analyzed the two-dimensional seepage in rubble dams using a one-dimensional non-Darcy seepage equation. Engelhardt and Finsterle [12] carried out hydraulic-temperature tests on broken rock mass mixtures as filling materials for underground storage rooms for nuclear wastes and estimated the permeability, thermal conductivity, specific heat capacity, and other parameters of different broken rock mass mixtures by the reversed modeling method. Miao et al. [13] and Ma et al. [14, 15] examined the seepage characteristics of broken rock mass with different porosities using self-made testing devices on broken rock mass. Bai et al. [16] proposed a type of plug model when studying the water bursting mechanisms of collapse columns. In this model, the porosity of the broken rock mass was found to vary with mass exchange, and the changing rate of the mass varied according to the given laws. Yao et al. [17] took into account the dissolution of water on broken rock mass, along with the transport of water for fine particles in broken rock mass, and established a type of broken rock mass deformation-water seepage coupling dynamic model. In this model, the mass of the broken rock mass was found to vary with the water dissolution and transport, which gave rise to changes in the porosity and permeability. Furthermore, there have also been some research studies conducted regarding the flow laws of the solid-liquid phases using numerical calculations. For example, Presho and Galvis [18] proposed a type of multiscale and generalized finite element method for changed flux fields in order to research the two-phase flow in porous media. Barba et al. [19] and Monaghan [20] studied the complex flow of solid-liquid using smoothed particle hydrodynamics. Hai and Hong-yang [21] carried out studies on the basic laws of water-sand two-phase flow in broken rock mass utilizing lattice Boltzmann method and established a numerical calculation model of water-sand two-phase flow in broken rock mass.

The water-sand two-phase seepage problems in broken rock mass are very complex. For the convenience of research, this study took the broken rock mass as seepage material and sand and water as the seepage fluids and assumed that there were no significant changes in the structures of the broken rock mass during the water-sand seepage processes; the pores between the broken rock masses were occupied with water and sand; the volume fraction, the mass fraction, and volume concentration also varied with time. This study established a water-sand two-phase seepage dynamic model of the broken rock mass based on such hypothesis and proposed a corresponding dynamic numerical computation method. Then, using FORTRAN language, a numerical calculation research was carried out on the seepage law of the water-sand two-phase fluids in the broken rock mass, which laid a theoretical foundation for the further study of the water-sand inrush mechanism during coal mining.

2. Water-Sand Two-Phase Flow Dynamic Model in a Broken Rock Mass

2.1. Motion Parameter Characteristics of the Water-Sand Two-Phase Flow

There are two phases of seepage fluids included in the water and sand phase of water-sand two-phase flow in broken rock mass. The space coordinates of any mass point in a broken rock mass will be identically equal to the material coordinates under the condition of not considering the structural deformation of the broken rock mass. Meanwhile, the space coordinates of the water-sand particles will vary with time. It was believed that any point in the space was occupied by both seepage material and fluid in poromechanics. In this way, there were three types of particles at any spatial point, which included broken rock mass particles, water particles, and sand particles, respectively. Their material coordinates were recorded as , , , , , , and , , , respectively. It was evident that the space coordinates , , and of broken rock mass particles were identically equal to their material coordinates , , and , while the space coordinates , , and of the water particles and space coordinates , , and of the sand particles vary with time.

In this study, since the deformation of the broken rock mass was not considered, its position and shape always overlapped with those of the control volume. The water and the sand were constantly flowing into the pores of the broken rock mass. Therefore, they had no certain overall position and shape and only partial position and shape, which corresponded to an infinitesimal element. Figure 2 illustrates the position and shape diagram of infinitesimal element, including the three types of infinitesimal media: broken rock mass, water, and sand at moment and moment . Figures 2(a) and 2(b) show the broken rock mass infinitesimal elements of the two moments, which had the same boundary and set of particles. However, the set of particles of the water and sand in its pores are obviously different. Therefore, strictly speaking, Figures 2(a) and 2(b) demonstrate the infinitesimal volume rather than the infinitesimal material. Figures 2(c) and 2(d) are the water infinitesimal element in moment and moment , respectively. In Figures 2(c) and 2(d), it can be seen that the set of particles of the water are the same in infinitesimal element. However, the spatial position was different. Figures 2(e) and 2(f) are the sand infinitesimal element in moment and moment , respectively. The set of particles of the sand were the same in infinitesimal element, but the spatial position was different.

The spatial positions of the broken rock mass infinitesimal element , water infinitesimal element , and sand infinitesimal element overlapped in moment . The position and volume of the broken rock mass infinitesimal element remain unchanged after experiencing time . The water infinitesimal element occupied spatial area , and its size and shape were different from . The sand infinitesimal element occupied spatial area , and its size and shape were different from .

Due to the fact that the movement of boundary was caused by the water flow, the material derivative of was determined to be related to the speed of the water particles :

Similarly, the material derivative of was related to the speed of the sand particles :

Given that the porosity between the broken rock masses was , the mass densities of the water and sand were and , respectively; and the volume fractions of the water and sand phase were and , respectively. Then, the mass concentrations of the water and sand were as follows:

The seepage velocities of the water and the sand were and , respectively, while the average velocities of the particles were and , respectively. Since the pores in the broken rock mass contained water and sand, the volume fraction was considered in the Dupuit-Forchheimer relation as follows:

2.2. Mass Conservation Equations

During the seepage process, the water-sand flow displayed small changes in mass density, which were considered and . Then, the water-sand phase mass conservation equations were, respectively, as follows:

2.3. Momentum Conservation Equations

In this study, under the condition of ignoring the capillary force, the momentum exchange of the water-sand was reflected in the relationship between the permeability parameter and volume fraction. In other words, the water phase permeability , non-Darcy flow factor , and acceleration coefficient were the functions of the volume fraction , and the sand phase permeability , non-Darcy flow factor , and acceleration coefficient were the functions of the volume fraction .

Water is considered to be a Newton fluid, and the momentum conservation equation is written on the basis of the Forchheimer relation as follows:where is the fluid pressure; is the dynamic viscosity of the water; the other parameters were the same as mentioned above.

Wet sand is considered to be a power-law fluid [22]. In this study, the consistency coefficient of wet sand was set as , the power exponent as , the apparent viscosity as , and the effective permeability as . Then, the momentum conservation equation was as follows:

Then, due to the relationship between the effective viscosity and apparent viscosity :

Equation (7) could then be written as

By setting the fluidities of the water and sand as and , respectively, the following could be obtained:

Then, according to (8) and (10), (7) and (9) could be simplified as follows:

2.4. Auxiliary Equations

Due to the fact that there were only two types of seepage fluids in the fractured rock pores, the water phase volume fraction and sand phase volume fraction met the following conditions:

The porosity evolution was determined by the porosity compressibility equation as follows:where represents the porosity when the pressure was ; is the pore compressibility coefficient; denotes the fluid pressure; represents the initial fluid pressure.

2.5. Permeability Parameters in the Water-Sand Two-Phase Flow in the Broken Rock Mass

During the water-sand two-phase flow in broken rock mass, the volume fractions of the water and the sand in the pores change with time. Therefore, fluidity, non-Darcy flow factor, and acceleration coefficient also change with time, and the change rule is affected by the sand grain size, graduation, and porosity. Therefore, the determinations of fluidity, non-Darcy flow factor, and acceleration coefficient are the key core elements to establishing the water-sand two-phase flow dynamic model in a broken rock mass.

According to the results of this study’s laboratory testing [22], it was found that the fluidity , non-Darcy flow factor , and the acceleration coefficient of water phase and the effective fluidity , non-Darcy flow factor , and acceleration coefficient of sand phase all conformed with the following polynomial relations:

From (14) to (19), represents the normalized volume fraction:where and are the left and right limits of the sand phase volume fraction, respectively; is the coefficient; , , and change with the porosity; ; the permeation parameters when and ( and ) could be calculated according to and , respectively.

The relationship between the equations and physical quantities in the above dynamic model is shown in Figure 3. As can be seen in the figure, link ① is the calculation of the effective viscosities and and the permeability parameters of the water and sand phases, which were mainly obtained using calculations (16) to (20) on the basis of the tests [22]. Link ② is the calculation of the volume fractions and according to the mass conservation equations (5) and the auxiliary equation (12). Link ③ is the calculation of the pore pressure according to the mass conservation equations (5) and the auxiliary equation (13). Link ④ is the calculation of the seepage velocities and according to momentum conservation equations (11).

3. Water-Sand Two-Phase Flow Dynamical Numerical Calculation in the Broken Rock Mass

3.1. Numerical Calculation Process

In this study, in order to further examine the water-sand migration rule in a broken rock mass, and in accordance with the abovementioned water-sand two-phase flow dynamic model, a numerical calculating program of the water-sand two-phase flow in a broken rock mass was programed using Fortran language program. Figure 4 illustrates its numerical calculating procedure. From the figure we can see that , , , , , , and were initially solved. This step was obtained through the basic calculations of the tests. At the same time, the effects of such factors as the grain sizes of the broken rock mass and the grain sizes of the sand should be considered. Since this step only involved algebraic operations, there was no need to construct algorithms. At this point, the numerical solutions of the sand volume fraction and the pore pressure were, respectively, obtained by the numerical integration of the mass conservation equation. The algorithm was roughly the same. Finally, the numerical solutions of the seepage velocities and were obtained through the numerical integration of the momentum conservation equation. For the occurrences of water-sand inrush disasters during coal mining processes, the focus is usually on the variation rule of the sand mass loss amounts Q under different influencing factors. In this study, the sand mass loss amounts could be solved in accordance with the variation law of the seepage velocity of the sand in the model. As the sand flowed out of the broken rock mass, the mass of the sand in the broken rock mass continuously decreased, and the volume fraction decreased, accordingly. During the time period of , the sand mass loss amounts flowing out of the broken rock mass werewhere is the sand mass density; is the seepage velocity of the sand; and denotes the width of the numerical model.

Through this study’s numerical calculation results, not only could the loss law of sand in the broken rock mass under the different influencing factors be obtained, the variation law of the permeability characteristics in the broken rock mass could also be obtained. Therefore, in this study, a theoretical foundation for the study of the mechanism of water-sand inrush disasters in shallow coal mining seams could be successfully laid.

3.2. Numerical Calculation Model
3.2.1. Establishment of the Model

In order to facilitate the study of the relationship between the sand grain flow, sand loss amounts and broken rock grain size, sand grain size, and porosity, it was proposed in this study to establish a cylindrical calculation model with a radius of  m and height of  m. The top interface of the model was used as the origin, and a plumb in a downward direction was used as the forward direction of Ox axis, in order to establish a numerical calculation coordinate system. Then, for the purpose of monitoring the change law of the sand grain flow in the model, the model was evenly divided into ten units, with a total of eleven nodes. The node spacing was 0.03 m, as shown in Figure 5. It should be noted that, in order to compare with the results of the test, in this numerical calculation, the sand grain and broken rock mass were evenly mixed. Then, water was continuously added from the top of the model, which formed the water-sand two-phase fluid. Also, during the calculation process, the sand was no longer added.

In order to prove the theoretical model’s reasonability with the experimental results shown in [22], the basic parameters of the model were set as follows: the mass density of the water was  kg/m3; dynamic viscosity of water was  Pa·s; mass density of sand grain was  kg/m3; consistency coefficient was ; power exponent was ; compressibility coefficient of the broken rock was  Pa−1; and porosity was under the reference pressure .

3.2.2. Numerical Calculation Scheme

(1) The variation law of the permeability parameters in the water-sand two-phase flow in broken rock mass has mainly been studied from the aspect of the variation law of sand phase seepage velocity and sand phase volume fraction , as well as the variation of such parameters as the sand flow phase , acceleration coefficient , and the non-Darcy flow factors during the flow process of the water-sand two-phase flow in broken rock mass. In this scheme, the particle sizes of the broken rock mass were between 12 and 15 mm; the sizes of the sand grains were from 0.074 to 0.25 mm; the top pressure was  Pa; the bottom pressure was  Pa.

(2) The factors, such as the particle size of the broken rock mass, sand grain size, and porosity, have important influences on inducing water-sand inrush disaster in a mine. Therefore, it has become very important to study the relationship between these factors and the amount of sand loss. This study’s detailed scheme was as follows: The influence of the particle size of the broken rock mass on the amount of sand loss: the sand sizes were between 0.074 and 0.25 mm; with a top pressure of  Pa, a bottom pressure of  Pa, and a porosity of . The particle sizes of the broken rock were 5 to 8 mm; 10 to 12 mm; and 12 to 15 mm. The influence of the grain sizes of the sand on the amount of sand loss: top pressure was  Pa; bottom pressure was  Pa, with a porosity of ; the particle sizes of the broken rock were 12 to 15 mm; the sand sizes were from 0.074 to 0.25 mm (fine sand); 0.25 to 0.59 mm (medium sand); and 0.59 to 0.83 mm (coarse sand). The influence of the porosity on the amount of sand loss: top pressure was  Pa; bottom pressure was  Pa; the particle sizes of the broken rock were between 12 and 15 mm. The sand sizes were from 0.074 to 0.25 mm, and the porosity was 0.565, 0.575, and 0.585.

4. Results and Discussion

4.1. Variation Law of the Permeability Parameters in the Water-Sand Two-Phase Flow in the Broken Rock Mass
4.1.1. Variation Law of the Sand Flow

The variation laws of seepage velocity and volume fraction of the sand phase are shown in Figure 6. The abscissa is the Ox axis of the model from top to bottom.

It can be seen from the figure that at  s, the seepage velocity of the sand was 0 m/s, which indicated that the sand had not been driven by water at the initial time. When  s and  s, the seepage velocity of the sand increased from the top of the model (0 to 0.03 m) to the lower part of the model (0.27 to 0.3 m) as a whole. However, the seepage velocity of the sand in the middle section of the model (0.06 to 0.24 m) remained essentially unchanged at  m/s, as shown in Figure 6(a). It can be seen that the phase volume fraction of the sand was evenly distributed in the model at  s, which accounted for 40% of the total volume. With the increasing seepage time, the volume fractions at the top of the model (0 to 0.03 m) and the bottom (0.27 to 0.3 m) were reduced to 0%, while the sand volume fraction also decreased with the increase of time in the middle of the model between 0.06 and 0.12 m. The sand volume fraction remained basically unchanged between 0.12 and 0.24 m, which was determined to be not related to the seepage time, as shown in Figure 6(b).

4.1.2. Variation Law of the Permeability Parameters during the Sand Migration

As can be seen from the figure, the permeability parameters also showed large differences during the migration of the sand grains along the Ox axis. The sand fluidity actually reflected the change of the sand phase permeability during the sand migration process. At the top of the model, the sand grains were moving into the middle section, due to being driven by water. The sand content at the top became increasingly less, which resulted in a relative increase in the sand phase permeability. Likewise, the sand grains at the bottom of the model also flowed out quickly, and sand phase permeability at the bottom increased. However, due to the sand in the middle section being supplemented by the sand at the top, the result was that the sand phase permeability was almost unchanged in 60 s. Then, with the increase in time, the increasing range of the permeability continued to expand into the center of the model, as shown in Figure 7(a). The acceleration coefficient characterized the inertia force during the sand migration process. It was observed that as the volume fractions of the sand at the top and bottom decreased, the quality also decreased. The acceleration coefficient also showed the permeability characteristic of first being small at the top and bottom and basically unchanged in the middle. The acceleration coefficient of some nodes in the central part also displayed a decreasing trend with the increase in time, as shown in Figure 7(b). The non-Darcy flow factor reflected the nonlinear feature of the sand’s migration process, in which the higher the numerical value was, the more obvious the nonlinear characteristic would be. The degree of sand nonlinearity was found to be highest at the initial time, and the nonlinearity decreases rapidly with the loss of sand at the top and bottom. Meanwhile, in the middle section, the change of the volume fraction of the sand was not obvious, while the non-Darcy flow factor displayed little change and also had high nonlinear characteristics, as can be seen in Figure 7(c).

4.2. Analysis of the Main Influencing Factors of the Water-Sand Inrush Disasters
4.2.1. Effects of the Particle Sizes of the Broken Rock Mass

Figure 8(a) shows the numerical calculation results of the variation of the sand loss amounts with different grain sizes in the broken rock mass. The following can be seen in the figure. (1) As the grain sizes of the broken rock mass increased, the sand loss amounts also increased, and the mine became more prone to water-sand inrush disasters. (2) The relationship between the grain sizes of the broken rock mass and the sand loss amounts was nonlinear. When the particle sizes of the broken rock mass ranged from 12 to 15 mm, the sand loss amounts were 882 g in 60 s. When the particle sizes of the broken rock mass ranged from 10 to 12 mm, the sand loss amounts were 624 g. When the particle sizes of the broken rock mass ranged from 5 to 8 mm, the sand loss amounts were 487 g. When the particle sizes of the broken rock mass were between 12 and 15 mm, the increases of sand loss amounts were larger than those of the previous two particle sizes. Figure 8(b) displays the results of this study’s indoor experiment [22]. It can be seen in the figure that as the particle sizes of the broken rock mass increased, the sand loss amounts also increased, which was the same as the numerical calculation. However, the total sand loss amounts were determined to be smaller than the numerical calculation results. For example, when the particle sizes of the broken rock mass were between 12 and 15 mm, the total sand loss amounts were 489.9 g. When the particle sizes of the broken rock mass were between 10 and 12 mm, the total sand loss amounts were 298 g. When the particle sizes of the broken rock mass were between 5 to 8 mm, the total sand loss amounts were 15 g. This was mainly due to the fact that the test system needed to achieve the set porosity by displacement loading, and the gravel grains would be extruded and deformed during the process of the displacement loading. This resulted in major changes in the gradation of the gravel, and the original pores became filled with the squeezed crushed grains. This caused the actual pores to be less than the theoretical calculation value, and the numerical calculation had effectively avoided this error.

4.2.2. Effects of the Sand Sizes

Figure 9 shows this study’s numerical calculation results of the variations of the sand loss amounts for the different sand sizes. The following can be seen in the figure. (1) During the entire seepage process, the speed of the sand losses during the first 10 s was fast, and the losses of the three sand sizes displayed only slight differences (259 g, 251 g, and 231 g, respectively). Also, the differences among sand loss amounts with the three sand sizes became wider with the continuous increases in the seepage time. (2) The sand with grain sizes from 0.074 to 0.25 mm had total loss amounts of 437 g in 60 s. The sand with grain sizes between 0.25 and 0.59 mm had total amounts of 398 g in 60 s. The sand with grain sizes between 0.59 and 0.83 mm had total amounts of 333 g in 60 s. Under the same condition, with the increases in the sand sizes, it was found that the smaller the sand loss amounts were, the easier it became for the water-sand inrush disasters of fine sand driven by water to occur in the mining of shallow coal seams. Therefore, prevention methods should be emphasized in the mining processes, and the drainage of loose sand grain aquifers should be made in advance, in order to prevent fine sand being acted upon by flowing water.

4.2.3. Effects of Porosity

It can be seen in Figure 10 that the sand loss amounts with different porosities had the following rules:

(1) The speed of the sand grain loss in 10 s was significantly higher than that in the later periods. That is to say, the slope was , which was mainly due to the fact that the sands at the bottom of the calculation model were close to the free face, and the sands, under the action of the water, flowed out very quickly and easily. As the time increased, the sand grains in the upper part of the model flowed into the middle section. In the cases of constant porosity, the sand grains in the middle were blocked, which led to slow losses of the sand grains. At the same time, it was found that there was a matching relationship between the sand sizes and the porosity. By taking as an example, the sand loss amounts became reduced after 10 s. However, it was still higher than that of the other sand loss amounts with two low porosities.

(2) With the increases in the porosity, the sand loss amounts became larger. For example, in the case of , the total amount of sand loss was 889 g in 60 s. Meanwhile, at and , the total amounts of sand losses were 455 g in 60 s. These results indicated that the larger the porosity of the broken rock mass was, the greater the probability of the water-sand inrush disasters and the amount of sand inrush would be. These findings were determined to be similar to the rule in the experiment, as well as the conclusions of the field observations [22].

5. Conclusions

In this research study, a basic examination of the nonlinear behavior of water-sand two-phase flow in the broken rock mass was carried out using both theoretical analysis and numerical calculation. This study mainly included the effects of porosity, particle sizes of the broken rock mass, sand sizes, sand content, and the lithology on the permeability parameters. The main conclusions which were obtained through the results of this study’s research experiments were as follows:

(1) In accordance with the mass conservation law and the momentum conservation law of the water and sand phases, and the porosity compressibility equation of the broken rocks, basic equation, and partial auxiliary equations of the system were derived. A water-sand two-phase flow dynamic model in a broken rock mass was established, and the initial and boundary conditions of the model were discussed.

(2) An algorithm of dynamic response of the water-sand two-phase flow system in a broken rock mass was constructed. The relationships between the porosities, volume fractions, pore pressures, seepage velocities, and permeability parameters were examined. Furthermore, the influences of the particle sizes of the broken rock mass and the sand sizes on the permeability parameters were determined.

(3) A numerical calculation model of the water-sand two-phase flow in broken rock mass was established. By using this calculation model, the variation laws of the sand grain parameters (e.g., the seepage velocity , volume fraction ), as well as the variation laws of the permeability parameters of sand phase (the non-Darcy flow factor ; acceleration coefficient , fluidity , etc.), were obtained. Due to the fact that the sand grains were not replenished in the numerical calculation, the sand grain in the upper part of the model gradually decreased with the increases in the seepage time. Also, the lower part of the model became rapidly reduced due to its closeness to the entrance. However, the sand grains in the middle section of the model were observed to be relatively unchanged. At the same time, the corresponding law of the sand-phase permeability parameters was also produced. In other words, with the increases in the seepage time , the sand-phase permeability parameters gradually decreased in the upper and lower ends of the model and basically remained constant in the middle. Also, the sand-phase nonlinear characteristic was higher in the middle than in the end sections.

(4) The particle sizes of the broken rock mass, sand sizes, porosity, and other factors were found to have important influences on the water-sand inrush disasters. Also, the sand loss amounts were greatly increased with the increases in the particle sizes of the broken rock mass and porosity and the decreases in sand sizes.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This paper was supported by the Research Fund of the State Key Laboratory of Coal Resources and Safe Mining, CUMT (13KF02); the National Natural Science Foundation of China (51704095 and 51774110); the Research Fund of the State Key Laboratory of Water Resources Protection and Utilization in Coal Mining, Shenhua Group (SHJT-16-30.11); the Research Fund of the State Key Laboratory for Geomechanics and Deep Underground Engineering, CUMT (SKLGDUEK1310); the Research Fund of the Henan Key Laboratory for Green and Efficient Mining and Comprehensive Utilization of Mineral Resources (Henan Polytechnic University) (S201619); the Key Project of Science and Technology of Department of Education of Henan Province (no. 14A440001).