Physical Simulation of the Water-Conducting Fracture Zone of Weak Roofs in Shallow Seam Mining Based on a Self-Designed Hydromechanical Coupling Experiment System

Due to inappropriate mining practices, water-conducting fracture zones can develop in an aquifer, not only destroying the surfacewater environment but also causing water inrush, even hurting or killing workers. To avoid such disasters, investigating and simulating the evolution mechanism of water-conducting fractures are becoming a research focus in mining engineering, especially regarding the organisation and development of fractures. Our work mainly involved the design of low-strength analogous materials and the simulation of fracture evolution for weak-roof problems in shallow seam mining based on a selfbuilt experimental hydromechanical coupling system. The experimental results show that the vertical stress in the roof increases first as the working face approaches and finally decreases to near its initial value as the working face passes. The relationship between fracture depth and coal-seam excavation distance is obviously nonlinear. The leakage velocity of surface water remains stable in the early stage of excavation and increases when the fracture develops through the main aquifuge. The maximum fracture depth is 76.18m for the Yili coal mine with weak roofs and shallow coal seams. In addition, we numerically simulated and verified the evolution patterns with the FLAC3D platform. The simulated fracture depth of the Yili coal mine agreed with the in situ borehole observation very well and was more accurate than the output of the empirical formula. Our work provides new methods and relevant data for research on the evolution of water-conducting fractures in weak roofs during shallow seam mining.


Introduction
The development of water-conducting fracture zones in overlying strata is a common mining problem. Once a water-conducting fracture zone develops in an aquifer, groundwater quickly begins to leak, causing problems in the surface environment and serious economic loss [1]. Therefore, water movement analysis is essential to research on the mechanism behind water-conducting fracture zone development in overburden strata during mining. Many scholars have predicted the development of water-conducting fracture zones in coal mining [2][3][4][5][6][7][8][9][10]. Zhang [11] introduced a mechanism for water-conducting zone development and proposed improved measures in coal mining based on hydrogeological conditions. A numerical simulation was carried out and compared with a field test to evaluate the risk of water inrush through the roof in Wuyang coal mine [12]. The impact factors for the development of a water-conducting fracture zone in a water-rich roof produced by long-span mining were analysed by Bai and Tu [13]. Majdi et al. [14] predicted the evolution mechanism of the height of the distressed zone of a roof in long-wall mining and proposed a simple method to estimate the relationship between the caving zone of the roof and the thickness of the coal seam. The overburden model can be regarded as a beam (2D) or a slab (3D), and some scholars have applied a numerical simulation method to analyse the bending problem of this model [15,16]. Kong et al. [17] discussed the impact factors for the stress intensity and bearing capacity of roofs in long-wall mining based on theoretical analysis, numerical simulation and field tests. Wang et al. [18] studied the dynamic response mechanism of a mining roof with PFC2D.
A similarity model test is an effective method for researching the spatial relationship between geological and engineering structures [19]. Laboratory physical simulation of large-scale engineering problems based on similarity model tests is widely used in coal mining [20][21][22]. Chen et al. [23] developed similar low-intensity materials and applied them as mining backfill. Guo et al. [24] simulated water inrush from a mining floor and proposed preventive measures for coal mining.
International scholars have performed many studies on the stress and deformation of mining roofs. However, there are few studies on shallow seams with weak roofs and fragile surface environments, such as that in Xinjiang and other locations in Northwest China. When water inrush occurs in mining, surface water can be depleted, and a fragile ecological environment is vulnerable to damage [25]. Therefore, the study of fracture evolution and water inrush from weak roofs in shallow seam mining has great significance for resource development and ecologically sustainable construction.
Based on the geological conditions of shallow seams in the Yili coal mine, Xinjiang, a series of analogous materials with low strength and low permeability are developed. On this basis, the response of overburden strata and the evolution of water-conducting fractures in mining are simulated by a self-designed experimental hydromechanical coupling system. Moreover, the evolution mechanism of waterconducting fractures in weak roofs during shallow seam mining is analysed.

Conditions of Mine Wells and Work Faces
The mining field is located in the southeastern part of Huocheng County, Yili Kazakh Autonomous Prefecture, Xinjiang Uygur Autonomous Region. The mining field spans 18.74 km in the E-W direction and 9.32 km in the N-S direction, with an area of 113.3 km 2 . The project unit is the Yili New Coal Mine Co., Ltd., with a construction scale of 6 million tons per year. The 21103 working face is the first mining face, the floor elevation ranges from +686 m to +886 m, the ground elevation ranges from +853 m to +950 m, and the fully mechanised caving face is located in the south-eastern part of the area at +660 m. The design working length is 1150 m, the inclination length is 120 m, and the coal thickness is 3.8-4.2 m. The main seam of the 21103 working face is 21-1 coal with an average depth of 120 m. The dip in the stratum is between 4°and 8°, and the average dip is 6°. A comprehensive mechanised coal mining method with longwall retreat and a double-cylinder shearer is adopted to cut coal, and the full-up fall caving method is used to manage the roof. The distribution from the coal seam to Quaternary strata (selected) is shown in Figure 1.
A comprehensive borehole histogram is shown in Figure 2. An aquifer with weak pores and a medium water content in the Quaternary layer is the indirect source of water filling the working face through the fracture zone. The Paleogene aquifer with a water-rich characteristic is the direct source of water filling of the working face through the caving zone. The fault structure is not developed in the original overburden strata, the groundwater recharge is poor, and the hydrogeological boundary is simple.
Taken together, the geological structure of the mining area, borehole histogram, and in situ protolith test (Table 1) indicate a shallow coal seam with weak-thin overburden strata. Although the structure of the stratum in the mining area is relatively complete, a fracture zone can easily develop in the overlying strata due to improper mining methods. The surface water and shallow groundwater can flow through the fracture zone of overburden strata to the working face in mining, leading to disastrous water inrush and serious damage to the fragile surface environment due to poor groundwater recharge conditions.

Model Establishment and Experimental
Testing of Water-Conducting Fracture Zone Evolution  2 Geofluids load must follow certain rules [26,27]. The ratio of a physical quantity with the same dimensions in the prototype and the model is called the similarity scale, which is expressed as C. L is defined as length, γ is severe, δ is dis-placement, E is elastic modulus, σ is stress, ε is strain, σ t is tensile strength, σ c is compressive strength, c is cohesion, φ is friction angle, μ is Poisson's ratio, and f is friction coefficient.  3 Geofluids (1) All dimensions of 1 are equal to 1, and the same scale is equal to the scale of the same dimension.
(2) According to the equation for the constitutive equation, the relationship between the displacement similarity scale, the geometric similarity scale, and the strain similarity scale is satisfied.
(3) According to the dimensional analysis, the relationships among the stress similarity ratio C σ , bulk density similarity ratio C γ , strain similarity ratio C ε , and geometric similarity ratio C L are as follows: 3.2. Similarity Test Materials and Method. Regarding the geological conditions of the Yili coal mine, the gravel layer and surface silt layer are aquifers with numerous fractures, large porosity. and low strength and contain surface water and shallow groundwater. The third layer (mudstone) serves as a water barrier, while the fourth layer (coarse sandstone) and the fifth layer (muddy siltstone) act as supports for the overlying strata. On the basis of these condition and research on similar materials in fluid-solid coupling [28], similar materials with low strength and low permeability are designed. The main components of the aggregate are sifted fine sand (120 mesh) and barite. The binder is gypsum powder and Vaseline. The permeability is determined by the silicone oil content. To explore the impact of each component on the strength and permeability of materials, an orthogonal experiment with 5 factors and 5 levels was adopted (  (Table S-2). According to previous studies [29], the designed analogous materials meet the requirements of low strength 4 Geofluids and low permeability. Consequently, the proportions for similar geological models can be determined from the information in Table 2.

The Model Test System.
To effectively explore the mechanism of the development of a water-conducting fracture zone in overlying strata during shallow seam mining, we designed a two-dimensional experimental "hydromechanical" coupling system ( Figure 4). This system can simulate the evolution of water-conducting fractures with surface water and shallow groundwater in mining. Compared with the traditional similarity model experiment, this system can simulate more realistic conditions. According to the borehole histogram ( Figure 2), the original rock test data (Table 1), and the similarity theory, we set the geometric similarity ratio at 1 : 200 and the density similarity ratio at 1 : 1.5. In addition, the time similarity ratio is 1 : 20, and the stress similarity ratio is 1 : 300. The similarity model experimental device includes a test fixture, water-loading system, and measurement system. The length of the similarity model is 150 cm, the height is 64 cm, and the width is 20 cm. A geoconstruction similarity model is made layer by layer through material mixing, mould loading, compaction, air drying, etc. To reduce the influence of interlayer forces on structural stability, mica sheets are laid between layers. Ten stress sensors and three displacement lines ( Figure 5) are arranged to monitor the dynamic response of the stress in the overlying strata during mining. The stress measurement points are labelled No. 1 to No. 10 and are evenly distributed in the pillars in two rows, one 15 cm and one 30 cm above the coal seam. The main stress measurement systems are stress sensors and TST3826E static strain gauges. The three measurement lines are established at 7 cm, 27 cm, and 47 cm above the coal seam, denoted line A, line B, and line C, respectively. Displacement is mainly monitored by the displacement measurement line. Finally, to retain water without changing the permeability, a hydraulic-loading device is embedded into the top silt, and the edge is sealed with water-proof glue ( Figure 6). Because the rate of excavation of the original working face is 10 m per day (24 hours), the rate of excavation of the model working face, which can be calculated by the time similarity ratio, is 5 cm per 1.2 hours. Therefore, we excavate 10 cm per 2.4 hours to ensure that the overburden deformation is complete. Moreover, 15 cm protective coal pillars are set on both sides of the model to eliminate the impact of boundary perturbation. Therefore, the effective excavation distance of the model is 120 cm.

Testing Results and Data Analysis.
In the early stage of mining (Figure 7(a)), the rock pressure increases steadily. When the excavated length reaches 20 cm (40 m), the maximum displacement of line C is 0.05 cm (0.1 m), which is larger than those of lines A and B. The overburden structure is relatively stable. A weak separation layer appears in the direct roof. The development of the fracture zone is slow, and the height of the fracture zone is 0.5 cm (1 m). Meanwhile, the water level is stable (Figure 8).
At an excavation depth of 30 cm (60 m), the rock pressure continues to increase steadily. Fracture development gradually accelerates. The fracture zone has not yet developed to the aquifer. The water level does not decrease significantly. When mining proceeds to 50 cm (100 m) (Figure 7(b)), rock pressure monitoring shows that the stress at the No. 2 stress sensor decreases from the maximum value of 2.26 MPa to 1.15 MPa, as well as No. 6 stress sensor (Figures 8(a) and  8(b)). At the position of the open-off cut, the stress reaches the maximum value of 5.76 MPa before it decreases, as detected by the No. 9 stress sensor (Figure 8(c)). The displacement monitoring line shows that the maximum bending subsidence of line C is 0.82 cm (1.64 m). The height of the fracture zone has developed to 21.5 cm (43 m), and the fracture zone has reached the main aquifer and developed cross-layer fractures. Therefore, the water level drops rapidly beyond an excavation depth of 50 cm (100 m). With the development of the fracture zone, the growth rate of the fracture zone slows. In addition, a bent subsidence zone develops 46 m away from the coal seam. Due to the bent subsidence of the underlying strata, the overlying strata also possess a bearing capacity. Therefore, an obvious separation space occurs.
After excavation to 60 cm (120 m), stresses at the No.  (Figure 7(d)). The maximum height of the fracture zone develops to 34.65 cm (69.3 m), and the number of large cracks in the penetrating layer increases, which leads to the stratified conduction of the fracture zone. Meanwhile, a bent subsidence zone has developed and reached the aquifuge. A subsidence separation space is formed between the aquifer and the aquifuge. A large amount of surface water leaks into the separation space. Numerous tensile cracks occur in the aquifuge due to the large deformation. Thus, the aquifuge has lost its water-proofing ability, which results in the acceleration of seepage and a large amount of surface water loss.

Simulation of the Mechanism
According to the experimental results, the development of fracture zones is the most intuitive reason for water inrush. The crack propagates through the main aquifer, causing 5 Geofluids  The main mechanism is shown in Figure 9. The formation of a goaf causes the roof to lose support from the bottom. Under the influence of gravity and pressure from the overlying strata, the roof undergoes subsidence and bending deformation and the internal stress increases gradually. Stress concentration (Figure 9(a)) occurs at both ends and the middle of the roof above the goaf during mining. With the increase in stress concentration, tension and shear cracks occur because the stress in strata reaches the ultimate yield stress. In addition, the deformation of and stress concentra-tion in strata continue to increase. The cracks gradually propagate through the strata, resulting in the formation of fracture groups (Figure 9(b)). Then, a separation space is formed in the upper part of the curved strata, which causes the overlying strata results to exhibit the same mechanical behaviour (Figure 9(c)). Finally, a fractured zone gradually develops upward through the aquifuge to the aquifer, and then water inrush occurs. With excavation, caving gangue accumulates gradually and supports the overlying strata, resulting in decreased deformation of the overlying strata, and the stress in the strata cannot reach the ultimate yield stress. As a result, the existing cracks stop propagating, and  Figure 5: Distribution of stress and displacement monitoring points. 7 Geofluids no new cracks are generated. Therefore, the height of the fracture zone tends to stabilise instead of continuously increase (Figure 9(d)).

Theoretical Model and Equations.
To illustrate and verify the mechanism of water-conducting fracture development intuitively from a mechanical viewpoint, we simulate the excavation of the working face in the Yili coal mine by FLAC3D.
In FLAC3D, the Lagrange algorithm is used. From the Cauchy formula, the stress vector t on any inclined plane can be obtained as where σ ij is the symmetric stress tensor of a known point in the medium, representing the stress state of the point and n is the unit normal vector on any inclined plane.

Hydraulic load box Overburden similarity model
Overlying water Water-repellent coating

Geofluids
If the particle has a velocity v, then in an infinitesimal time dt, the medium will produce an infinitesimal strain determined by v i dt, and the corresponding strain rate component ξ ij is Its rotation rate component ω ij is where x is the vector representing the position of the particle. According to Newton's law of motion and Cauchy's principle of stress, if the particle acts on the stress σ ij and the physical force b i and has the velocity v i , there is Cauchy's equation of motion in a very small time span dt.
where ρ is the particle density. If the acceleration of the particle is 0, the above formula becomes the static equilibrium equation, Equations (6) and (7) are the geometric equation and equilibrium differential equation, respectively, so the constitutive equation describing the relationship between the stress and strain of the medium is also needed. The definitions in FLAC3D are as follows: where b σ ij is the rate of stress change, H is a known function, and κ is a parameter related to the loading history.

Geofluids
In FLAC3D, the constitutive model used in the elastic stage of the medium is a linear elastic model based on the generalised Hooke's law, including the isotropic elastic model, the orthotropic model, and the anisotropic model. In this paper, the isotropic elastic model is used in any single stratum, and its general expression is as follows: where E is the elastic modulus of the material, μ is Poisson's ratio, and G is the shear modulus. For the plastic stage, there are many constitutive models commonly used in geotechnical analysis, while for the excavation of rock and soil and the stability of the slope, the more commonly used model is the Mohr-Coulomb model [30]. Therefore, the Mohr-Coulomb model is selected as the plastic model of the numerical simulation in this paper. It mainly includes the shear failure criterion and tensile failure criterion; that is, where σ t is tensile strength, and N φ = 1 + sin φ/1 − sin φ.  Figure 11.

Model Establishment and Simulation
According to the plastic zone and stress distribution map, the plastic zone (Figure 11(a)) develops in the coarse sandstone layer at the top of the main aquifuge, with a height of approximately 78 m. The overburden structure experiences mainly tensile failure, which is reflected by the maximum principal stress (Figure 11(b)). The main reason for the failure of coarse sandstone layers is that the internal failure of tensile stress reaches the tensile strength limit of the layer, which is determined by the bending subsidence. The experimental results show that the generation and development of fractures on both sides and in the middle of the strata are consistent with the distribution of maximum principal stress (Figure 11(b)) and direction of shear stress (Figure 11(c)). With the excavation of the working face, the limit of the stress and plastic zone increases. The position of maximum vertical stress (Figure 11(d)) in overburden strata appears directly above the goaf. After the advancement of the working face, a pressure relief zone appears in the roof, and the    subsidence. As a result, the plastic zone penetrates the aquifer ( Figure 12). Overall, the simulation results are close to the experimental results, which verifies the correctness of the mechanistic analysis.

Correlation between the Excavation Distance of the Working Face and the Development of the Fracture Zone.
Controlling the maximum height of the water-conducting fracture zone is one of the main discriminant methods for water control in mining engineering. In engineering, to predict and prevent water inrush during mining, it is necessary to ensure that the maximum height of the waterconducting fracture zone does not exceed the height between the coal seam and the aquifer. At present, an empirical formula (Equation (12)) is usually used in engineering to predict the maximum height of the water-conducting fracture zone of weak overburden (compressive strength less than 20 MPa) in mining [31].
H f = 100M 5:1n + 5:2 where M is the thickness of the coal seam and n is the number of coal seams. For the working face in the Yili coal mine, Xinjiang, we apply the geological data obtained from field tests to empirical formulas. The thickness of the coal seam is 4 m, and the number of coal seams is 1. Therefore, the maximum height of the water-conducting fracture zone in long-wall mining is approximately 43.93 m.
Because of the inhomogeneity of the distribution of overlying strata, we also test the empirical formula (Equation (13)) for the maximum height of the water-conducting fracture zone in coal mines with medium-strength overburden (compressive strength between 20 MPa and 40 MPa).
The maximum height of the fractured zone is 61.43 m.
The experimental results reveal that the height of the fracture zone reaches 21.5 cm (43 m) when the working face advances to 50 cm (100 m). When the excavation distance reaches 70 cm (140 m), the height of the fracture zone is 34.65 cm (69.3 m), which is obviously higher than the height calculated by the empirical formula.
To quantitatively explore the relationship between the excavation distance and the height of the fracture zone, a nonlinear fitting is made between the excavation distance and the height of the fracture zone ( Figure 13). The results show that the fit converges where the Chi-square tolerance is 1e-9, R-square (COD) is 0.99975 and Adj. R-square is 0.99956, indicating that the fitting result is satisfactory.
From the fitting results, we find that the height of the water-conducting fracture zone in the overburden structure tends to be stable during mining, which is also consistent with the engineering results for another overburden in a shallow seam mine [32]. On the basis of calculations with increasing excavation distance, it is found that the height of the fracture zone steadily increases, and the maximum height of the water fracture zone is approximately 76.18 m.

Geofluids
On the basis of the observations from the borehole test (Table 3) of failure in the overlying strata in the 21-1 coal seam at the 21103-working face [33], the height of the caving zone is approximately 10.98-12.55 m, and the height of the fracture zone is 72.76-80.47 m. The average height of the fracture zone can be calculated to approximately 76.615 m. The field test results are very close to the experimental fitting results, but there is a large gap in the results of the empirical formula. While proving the validity of the experiment, the traditional empirical formula is not proven suitable for calculating the maximum height of the water-conducting fractured zone in the overburden in shallow seam mining.

Other Influencing Factors and Follow-Up
Study. According to previous studies, the factors impacting water inrush include the stress in the primary rock, the original overburden structure, and the water pressure [11]. This paper focuses on the unique Yili coal mine, Xinjiang, where the strength of overlying strata is weak and there is no extremely water-rich confined aquifer. Surface water and shallow groundwater are the main sources of recharge. Therefore, the different conditions of in situ stress, lithology, and confined aquifers that may exist in other shallow mining areas are not considered, and the study of other factors impacting the development of a water-conducting fracture zone in shallow seam mining and the universal laws governing this development are the focus of future work.

Conclusions
Physical similarity models are used to simulate the evolution mechanism of a water-conducting fracture zone according to the coal seam-overlying aquifer arrangement in the Yili coal mine, Xinjiang. The water-conducting fracture zone gradually expands with increasing mining distance. Moreover, the spreading speed tends to gradually increase, showing an obvious nonlinear distribution. The critical water-conducting fracture zone is reached at an excavation depth of 60 cm (120 m). Moreover, the depth of the fracture zone reaches 31 cm (62 m) and begins to penetrate the main aquifer. When the mining proceeds to 70 cm (140 m), the fracture completely penetrates the main aquifer of coarse sandstone, causing serious instability in the subsidence of the upper mudstone strata and a sharp decrease in the water level.
We also analyse the main mechanism of the evolution of the water-conducting fracture zone. After mining, the roof bends under the action of overlying strata and its own gravity. When the ultimate stress is exceeded, the fracture develops. We use FlAC3D to simulate coal seam excavation. The simulation results verify the correctness of the stress distribution in the testing and analysis. As mining continues, the fracture zone tends to stabilise and no longer grow primarily because the overlying strata are gradually supported by the collapsed strata. As a result, the overlying strata are no longer deformed, and the internal stress strength gradually fails to reach the ultimate stress strength.
By fitting the relationship between the depth of the fracture zone and the excavation distance, we find that the experimental results are similar to the field drilling test results. This proves the correctness and feasibility of the experimental methodology and demonstrates that the traditional empirical formula requires modification for calculation of the maximum depth of the water-conducting fracture zone in shallow seams in Northwest China.

Data Availability
In our manuscript, the experiment and test data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
No conflict of interest exists in the submission of this manuscript.