Impacts of longwall mining speeds on permeability of weakly cemented strata and subsurface watertable: a case study

Abstract The negative impacts of underground longwall coal mining on water resources and ecological environment are associated with changing overburden permeability and subsurface watertable. For weakly cemented strata characterized by short diagenesis, low strength, easy swelling and mudding with water, overburden permeability and subsurface watertable are highly sensitive to the advancing speed of the longwall. In this paper, the longwall panel 21103 in Yili No.4 Coal Mine, Xinjiang, China, where weakly cemented strata exist in Jurassic and Palaeogene systems, was taken as the case to study the response of overburden permeability and subsurface watertable to different advancing speeds. Results show that there is a butterfly-shaped permeability-increasing area comprised of a trapezoid and double V-shaped subareas in overburden. Generally, the permeability curves of aquiclude are of double peaks. In shallow aquifer, there is an exponential decline in the maximum drawdown of watertable, while a linear decline in lateral influence range as the speed goes up. If constrained by the ultimate watertable drawdown still being able to sustain the regular growth of vegetative covers, the optimal minimum advancing speed of panel 21103 is determined as 8 m/d, which is demonstrated by field verification.


Introduction
With national energy development strategically shifting to the west, a great deal of coal mines in Northwest China has become the mainstay in outputting coal resources. Locally, weakly cemented strata extensively occur in Jurassic and Palaeogene systems, and such young sedimentary strata are of short diagenesis, low strength, weak cementation, and easy swelling and mudding with water . There is no significant key stratum in weakly cemented strata and overburden tends to crack and collapse just following the advance of longwall coalface (Zhang, Sun, et al. 2019). Also, arid and semi-arid climate in Northwest China make local area scarce in freshwater and resistant ecological environment and accordingly, shallow water plays an irreplaceable role in sustaining ecological stability of such mining regions. In fact, inappropriate mining operations can lead to excessive development of water flowing fractures in overburden and loss of water-resisting capacity of aquiclude, which will be highly likely to trigger a series of problems like groundwater loss, land desertification and vegetation degradation sequentially (Rashid et al. 2014;Li et al. 2015). The contradiction between large-scale mining activities and ecological environment protection restricts the exploitation of coal resources over a long period. Existing studies show that advancing speed of coalface is considered a key factor that determines overburden deformation and movement (Liu 2014;Yang et al. 2015). In turn, such deformation and failure directly lead to the change of overburden permeability, and further cause the watertable decline. Apart from China, arid and semi-arid coal mining areas commonly exist in the United States, Australia, India, Pakistan, Iran et al., which makes the negative impact of mining activities on shallow water resource a worldwide problem especially in ecologically fragile area (Anawar 2013;Pazand and Sarvestani 2013;Krzyszowska Waitkus 2018;Gupta and Sharma 2019;Ayub and Ahmad 2020). Based on the lithological property of weakly cemented strata, the impacts of advancing speed on overburden permeability and subsurface watertable are analyzed in this study. This work is of great significance for underground coal mining in ecologically fragile areas both in China and all over the world.
Up to now, academics have escalated a good deal of researches on permeability variation of overlying strata subjected to mining disturbance, and it is commonly recognized that such a variation is directly related to strata movement, fracture development, stress state or even porosity change of overlying strata (Adhikary and Guo 2015;Meng et al. 2016;Khanal et al. 2019;Zhai et al. 2019). For instance, Wachtel et al. believed that coal mining triggered the change of overburden strain and further led to the changes of overburden permeability (Wachtel et al. 2014); Das thought that it was the upward development of mining-induced fractures that resulted in a significant increase of overburden permeability (Das 2000); Karacan and Wang et al. thought that, with coalface advancing, collapsed rocks in goaf gradually reconsolidated, which made for geo-stress redistribution of overlying strata and the change of overburden permeability thereof (Karacan et al. 2007;Wang et al. 2016); Poulsen, Nourani and Carcione et al. revealed that overburden permeability directly related to rock porosity, and a method to assess the change of overburden permeability was put forward based on the porosity parameters (Poulsen et al. 2018;Nourani et al. 2019;Carcione et al. 2020). Yang and Ma et al. found that there was a negative correlation between overburden permeability and surrounding rock pressure. And the overburden permeability can be effectively improved by borehole pressure relief (Yang et al. 2011;Ma et al. 2020). It can be seen that some basic laws in explaining overburden permeability change induced by mining operations and relevant evaluation methods have been proposed. However, there have been few studies associated with the impacts of advancing speeds on overburden permeability in weakly cemented strata, especially when lithological properties of such rocks are considered.
In addition, the impacts of mining operations on water resources attract lots of attention (Arkoc et al. 2016;Luan et al. 2020;Sahoo and Khaoash 2020;Ol ıas et al. 2021). For example, Tiway and Santana et al. pointed out that water resources were the most vulnerable environmental compartment directly affected by mining operations, especially in semi-arid mining areas (Tiwary 2001;Santana et al. 2020); Obiadi and Shepley et al. thought that mining-induced factures would serve as the main channels for water loss and pollutant inflow in aquifer (Shepley et al. 2008;Obiadi et al. 2016); Howladar investigated aquifer watertable and revealed that its change right correlated with the location and depth of coal mines, as well as the recharge and drainage conditions of aquifer (Howladar 2013); Sidle et al. thought that underground longwall mining could affect the length of cascades, extent glides, pool length, pool volumes and channel geometry of mountain stream channel (Sidle et al. 2000). Guo and Ping et al. pointed out that mining-induced fractures were the main factor causing the decrease of surface runoff (Guo et al. 2017;Ping et al. 2017). Existing studies mainly relate to the impacts of mining operations on water loss and water quality deterioration from perspectives of mininginduced fractures, mining depth, and recharge and drainage conditions of aquifer. In fact, for weakly cemented strata that attract great attention over recent years, there are few achievements about how advancing speed affect subsurface watertable based on the variation characteristics of overburden permeability.
In this paper, weakly cemented coal measure strata in Yili No.4 Coal Mine, Xinjiang, China are taken as the case to analyze the influence of advancing speed on overburden permeability and subsurface watertable. Yili No.4 Coal Mine is situated in arid and semi-arid areas. Because of insufficient precipitation, local ecological environment is relatively vulnerable to mining disturbance. As a matter of fact, the shallow aquifer plays an important role in maintaining the ecological stability of the mining area. Therefore, for longwall coal mining in similar areas like this, an appropriate advancing speed of coalface is necessary so that the watertable drawdown of shallow aquifer can be controlled within an allowable scope. Based on the strain-permeability function of weakly cemented rocks in Yili No.4 Coal Mine, the response of overburden permeability to different advancing speeds is analyzed. Then, a groundwater seepage model is developed and then used to investigate the impacts of different advancing speeds on subsurface watertable. The optimal advancing speed of panel 21103 in Yili No.4 Coal Mine is subsequently determined by taking the ultimate watertable drawdown that could still support the regular growth of vegetative covers as constraint condition. Finally, the rationality of optimal advancing speed of panel 21103 is verified by industrial tests on site.

Geological setting
Yili No.4 Coal Mine is located in Yining City, Xinjiang, China. The climate of the mining area is arid and semi-arid with annual precipitation ranging from 186.3 mm to 387.7 mm, and the ecological environment is recognized fragile. Coal measure strata mainly belong to Jurassic system and are of young sedimentary age, low strength, weak cementation and rocks thereof are naturally rich in clay minerals such as montmorillonite and kaolinite . The shallow aquifer is Quaternary gravel layer and just below it is the Palaeogene mudstone aquiclude. At the same time, the aquifer is mainly supplied by rainfall and snowmelt. Panel 21103 serves for pilot production as the first mining panel with length and mining height reaching 115 m and 3.5 m respectively. When panel 21103 is mined, solid coal is around the coalface. Here, the depth of coal seam is 55-150 m and 115 m on average. The dip angle of coal seam is 4-10 and 6 on average. The shallow gravel aquifer is covered by 20 m thick silt layer. Such an average 10 m thick aquifer plays an important role in maintaining the stability of ecological environmental in the mining area. In addition, the average thickness of mudstone aquiclude is 5 m. The average distance between coal seam and gravel aquifer is 85 m. The comprehensive geological column and geometry of panel 21103 are shown in Figure 1.

Mining situation
The initial advancing speed of panel 21103 is 4.8 m/d. There was a quite interesting phenomenon in the mining process of coalface. On-site observation shows that the water inflow from the gravel aquifer was low when such an advancing speed 4.8 m/d was adopted in consecutive coal mining, while the water inflow saw a sharp increase when the normal mining process was restarted from safety inspection or equipment maintenance. The bucket volume method is used to monitor the water inflow of coalface. During the advancing process, the water inflow monitoring points are arranged in the two retreat roadways adjacent to the panel 21103 coalface. And the water inflow of coalface is monitored by monitoring the water inflow in the bucket within the unit time. On-site observation data of water inflow changing with advancing speed are shown in Figure 2. It can be seen from Figure 2 that the average water inflow into panel 21103 was 11.7 m 3 /h when the advancing speed of the coalface was 4.8 m/d. But this water inflow went up to 83.0 m 3 /h (on December 8, 2016), 7 times higher than the routine level, after a 3-day production halt due to safety inspection (from December 3, 2016 to December 5, 2016). Such a high inflow event lasted for a week, causing water loss in quantity of shallow gravel aquifer. Similar phenomenon then occurred from January 2, 2017 to January 5, 2017. The maximum water inflow of coalface was 78.0 m 3 /h after a 2-day production halt due to safety inspection. When the advancing speed was reduced to 0.8 m/d due to the Chinese Lunar New Year holiday (from January 27, 2017 to February 6, 2017), the water inflow of coalface increased significantly. During the holidays, the maximum and average water inflow of coalface were 36 m 3 /h and 33.0 m 3 /h, respectively. Compared with the normal mining process, the average water inflow of coalface at 0.8 m/d advancing speed was 2.8 times of that at 4.8 m/d advancing speed. This shows that under the condition of weak cemented strata, the advancing speed has significant influence on the water inflow into coalface, and there is a negative correlation between them.

Impacts of advancing speeds on overburden permeability
In order to analyze the impacts of advancing speeds on overburden permeability, a solid-fluid coupling numerical calculation model, in accordance with the geological conditions of panel 21103, is established by means of FLAC 3D software (Zarlin et al. 2012;Do et al. 2015;Zhu et al. 2020). Based on the FLAC 3D model, the variation of overburden permeability under different advancing speeds is analyzed.

Model establishment
According to the comprehensive geological column (shown in Figure 1) of panel 21103 and the physical and mechanical parameters of coal and rock masses (listed in Table 1), the FLAC 3D numerical calculation model is established as shown in Figure  3. The model is 400 m, 315 m and 153.5 m in length, width and height respectively. In detail, the bottom and four lateral boundaries are fixed on displacement, while the top boundary is set as free. The length, width and height of panel 21103 are 200 m, 115 m and 3.5 m respectively, which are well in accordance with actual circumstances. The model is composed of cube grids with 5 m in both length and width. At the same time, the height-width ratio of the cube grids is 0.5-1.0 (shown in Figure 3). The total number of cube grids is 1.26 Â 10 6 . As mentioned above, panel 21103 is the first mining panel which is surrounded by coal body. In order to eliminate the influence of boundary effect of FLAC 3D model, 100 m wide protective coal pillars are set around panel 21103 both longitudinally and transversely. In order to make the numerical model converge quickly and reduce the computational effort, the coal-rock mass is assumed to be isotropic equivalent porous medium (Ebigbo et al. 2016;Salmi et al. 2017;Madadi et al. 2018;Zertsalov et al. 2020). The Mohr-Coulomb model is selected as the failure criterion of coal and rock masses and the seepage mode is set as isotropic. The solid-fluid coupling governing equations (Such as Transport Law, Balance Laws, Constitutive Laws and Compatibility Equations) used in the numerical model is the default equations in FLAC 3D software. See the manual of FLAC 3D software for details (Itasca Consulting Group, Inc. 2012).
In order to make the established FLAC 3D model reflect the field conditions of panel 21103 more accurately, the functions representing the relationship between volumetric strain and permeability of overlying strata are input into the numerical calculation model. y ¼ ð1:17 þ 8:23e À609:76x Þ Â 10 À19 0 x<2:6 Â 10 À3 y ¼ ð1:29À11:99e À952:38x Þ Â 10 À18 x ! 2:6 Â 10 À3 ðMudstoneÞ (1) where y is the permeability (m 2 ) and x is volumetric strain in both functions (1) and (2). In fact, the functions are derived from our previous study in which the permeability of weakly cemented mudstone and sandstone samples drilled from Panel 21103 were tested along a designed complete stress-strain path and resultant data about volumetric strain and permeability of weakly cemented rocks were obtained (Fan et al. 2020). In order to make the fitting functions more accurate, such data are fitted by virtue of piecewise function based on the original data and the results of volumetric strain-permeability function is shown in Figure 4. It can be seen from Figure 4 that, with volumetric strain going up, the permeability of both weakly cemented rocks decreases and then increases exponentially. R 2 , the goodness of fit, of each function is greater than 95%, indicating that the actual relationship between volumetric strain and permeability of weakly cemented rocks in complete stress-strain process can be well represented by the above functions.

Model calibration
In order to make sure whether the model is reasonable or not, a surface subsidence measuring line, highlighted in Figure 3, is arranged on the upper surface of the model along the layout direction of the panel 21103 coalface. After the model is meshed appropriately, the computational stage, where the coalface advances at the speed of 5 m/d, then undergoes and the resultant subsidence is output and exhibited in Figure  5, in combination with corresponding values measured on site. It can be seen from Figure 5 that both curves share a similar trend in general and there is also a good agreement between values themselves. Compared with the maximum subsidence (3.10 m) measured on site, such an index obtained by simulation is 3.03 m and the error between both values only accounts for 2.3%. As a result, it is believed that the physical and mechanical parameters of coal and rock masses (listed in Table 1) used in the FLAC 3D numerical model can genuinely represent the actual geology conditions of panel 21103.

Modelling scheme
Based on the solid-fluid coupling numerical simulation model established for panel 21103, five independent FLAC 3D numerical models are calculated at the advancing speeds of 3 m/d, 5 m/d, 7 m/d, 9 m/d and 11 m/d respectively, and then the variation of overburden permeability in weakly cemented strata under different advancing speeds is analyzed. In order to realize the difference in speed of coal mining, single excavation step spacing of each model is assigned as 3 m, 5 m, 7 m, 9 m and 11 m respectively. After each single excavation, the time from model calculation to automatic balance is counted as one day until all excavation is completed.

Modelling results
Vertical profiles of overburden permeability under different advancing speeds are shown in Figure 6. The acquisition time of permeability profiles is when the models are fully excavated and calculated to the automatic balance. What's more, the positions of permeability profiles are located in the middle of the models and are arranged along the panel coalface, which is exactly the same as the layout of the monitoring line of surface subsidence in Figure 3.
From Figure 6, it can be seen that the distribution characteristics of overburden permeability in weakly cemented strata are basically the same under different advancing speeds. Along the layout direction of the panel coalface, there is a butterflyshaped permeability-increasing area, which is composed of a trapezoid and double V-shaped subareas. Among them, the trapezoid permeability-increasing subarea is located in lower overburden closely above the goaf, and the double V-shaped permeability-increasing subareas are located on both sides above the trapezoid and symmetrical about the centreline of the panel. With the increase of advancing speed, the height of the trapezoid subarea linearly declines (shown in Figure 7) and such a change can be expressed by y ¼ -3.51x þ 60.22 with R 2 levelling at 0.98 (where y is referred to as the height of the trapezoid subarea, m; x referred to as the advancing speed, m/d; R 2 referred to as the goodness of fit). When panel coalface advances at the speed of 3 m/d, the height of trapezoid subarea peaks at the value of 49.5 m and by this time, its top has been at the shortest distance (30.5 m) from the bottom of aquiclude.
In addition, according to the simulation results shown in Figure 6, the lower the advancing speed, the worse the water-resisting ability of aquiclude. At the same time, the dominant seepage channels of aquiclude are located at the position corresponding to both ends of the panel coalface (shown in Figure 6a,b). When advancing speed is 11 m/d, the aquiclude is kept almost integral and just subjected to bending deformation (shown in Figure 6e). When advancing speed is 9 m/d, the effective water-resisting thickness at both sides of aquiclude gets thin (shown in Figure 6d). The effective water-resisting thickness would get thinner and the residual thickness being able to support effective water-resisting function is less than 50% of the original thickness of aquiclude when advancing speed is 7 m/d (shown in Figure 6c). Aquiclude would lose more water-resisting ability while penetrated by cracking regime when the advancing speed is 5 m/d (shown in Figure 6b). When the advancing speed is 3 m/d, the aquiclude is failure and the permeability increases sharply (shown in Figure 6a).
Aquiclude permeability at different advancing speeds is output and depicted as shown in Figure 8a. Permeability curves are of double peaks that are respectively located above both ends of panel coalface. With the increase of advancing speed, the maximum permeability of aquiclude decreases obeying an exponential function. It can be expressed as y ¼ (9.98e À0.24x þ0.45)Â10 À17 (where y is referred to as the peak of permeability, m 2 ; x is advancing speed, m/d), and the goodness of fit equals R 2 to

Impacts of advancing speeds on subsurface watertable
In order to analyze the impacts of different advancing speeds on subsurface watertable, a hydrological model for groundwater transport is established by GMS (Groundwater Modelling System) in accordance with the actual hydrogeological conditions of panel 21103. Based on the GMS model of panel 21103, the variation of subsurface watertable under different advancing speeds is analyzed.

Model establishment and scheme
This established GMS model for groundwater migration is of the same geometry as the FDM-based FLAC 3D model in Section 3.1, in which panel length and width are 200 m and 115 m respectively and at the same time, 100 m wide coal pillars around the panel are maintained. A cross section of the GMS model along horizontal gravel aquifer is shown in Figure 9. Four lateral sides of the model are assigned the same head as recharge boundaries and on such basis, absolute change of aquifer watertable in response to different advancing speeds is investigated without considering irregular recharge and drainage around the panel. Besides, the original water head of gravel aquifer is set as 10 m according to the actual hydrogeological conditions of panel 21103. Other strata below gravel aquifer are set to water saturated state before model calculation. What's more, the steady state flow model is selected for the GMS model.
As GMS is just a water flow software, it is cannot set different advancing speed in the model directly. However, for the established GMS model, strata permeability is a key parameter for groundwater migration calculation. The permeability of overburden under different advancing speeds can be obtained from Section 3.2. Assign the permeability parameters of overburden under different advancing speeds to five independent GMS models. As the permeability parameters themselves have time scales, the different advancing speeds (3 m/d, 5 m/d, 7 m/d, 9 m/d and 11 m/d) for GMS models are set indirectly by this method. In addition, the mined coal seam with an area of 115 m Â 220 m is set as pumping well so as to simulate aquifer water loss induced by underground coal mining. The drainage of the pumping wells under different advancing speeds is determined according to the actual water inflow of panel 21103 coalface at the corresponding advancing speed. Water inflow under different advancing speeds measured in panel 21103 coalface during normal continuous mining process is shown in Figure 10. Based on unary linear fitting, it is found that the relationship between actual water inflow and advancing speeds can be expressed as y ¼ -1.49x þ 18.98 (where y is referred to as actual water inflow in workface, m 3 /h; x is referred to as advancing speeds, m/d), whose goodness of fit R 2 is 0.75. According to this function relationship, drainage of pumping well under different advancing speeds is calculated and shown in Table 2.

Modelling results
The vertical drawdown of aquifer watertable is computed and exhibited in Figure 11.
It can be seen from Figure 11 that the vertical drawdown contour of gravel aquifer around panel is of oval profile and the maximum is located in the middle of the goaf. When panel 21103 coalface advances at the speed of 3 m/d, aquifer watertable decline reaches À4.83 m to the maximum extent, accounting for 48.3%. As advancing speed goes up, the maximum vertical drawdown shows significant decrease obeying an exponential function y ¼ -5.68e À0.21x -1.76 (where y is referred to as the maximum vertical drawdown of aquifer watertable, m; x is advancing speeds, m/d; the goodness of fit equals to 0.99). See Figure 12 for details. If taking À2 m vertical drawdown as a reference, the lateral influence range reaches 129.34 m, 109.08 m, 100.52 m, 77.7 m and 57.5 m when the advancing speed is 3 m/d, 5 m/d, 7 m/d, 9 m/d and 11 m/d respectively. Besides, the lateral influence range of À2 m vertical drawdown linearly shrinks with the increase of advancing speeds and the relationship between both variables obey a linear function y ¼ -8.75x þ 156.10 (where y is the lateral influence

Mechanism analysis and on-site verification
According to the above simulation results, the changes of overburden permeability and shallow aquifer watertable at different advancing speeds are obviously different. Overall, the higher the advancing speed, the smaller the increase of overburden permeability and the drawdown of aquifer watertable. In order to explain this interesting phenomenon, the characteristics of overburden caving and fracture development under different advancing speeds are analyzed, and then the influence mechanism of advancing speed on overburden permeability and aquifer watertable is revealed. On this basis, the advancing speed of panel 21103 in Yili No.4 Coal Mine is optimized, and the optimization result is verified by field test.

Mechanism
The essence influence of advancing speed on overburden permeability and aquifer watertable is the influence of advancing speed on overburden movement and deformation. A large number of studies have shown that after underground coal mining, overburden strata collapses and forms stable voussoir beam structure (Sofianos and Kapenis 1998; Tsesarsky 2012; Li et al. 2016). The voussoir beam is subjected to horizontal compressive force, which has a bearing capacity on the overlying strata. The voussoir beam structure and fracture development at different advancing speeds are shown in Figure 13. Research shows that the higher the advancing speed, the larger the break distance of single voussoir beam block (L 1 in Figure 13a), the bigger the horizontal compressive force (T 1 in Figure 13a), and the stronger the bearing capacity of voussoir beam structure (P 1 in Figure 13a) . But when the advancing speed is low, the results is just opposite. The lower the advancing speed, the shorter the break distance of single voussoir beam block (L 2 in Figure 13b), the smaller the horizontal compressive force (T 2 in Figure 13b), and the weaker the bearing capacity of voussoir beam structure (P 2 in Figure 13b). Strong supporting capacity of voussoir beam (P 1 in Figure 13a) at a high advancing speed can control overlying strata to produce slight mining fractures, and thus keep aquiclude in a good integrity (shown in Figure 13a). While weak supporting capacity of voussoir beam (P 2 in Figure 13b) at a low advancing speed will result in a large number of transfixion cracks in the overburden and poor integrity of the aquiclude (shown in Figure 13b). The higher the coalface advanced, the shorter time the cracks propagated, the faster the cracks closed, and the more gently the surface subsided (Zhang et al. 2009). In addition, for the weakly cemented strata, the rich clay minerals such as montmorillonite and kaolinite can be expanded and argillated with water, which can block some tiny mining cracks and has a good closing effect on the slight cracks produced at a high advancing speed . This explains why the higher the advancing speed, the lower the overall permeability of overburden and the lesser the aquifer watertable drops. Therefore, within a certain range, improving the advancing  speed of panel coalface can reduce the development of overlying strata fractures, reduce the increase of overlying strata permeability, reduce the decline of aquifer watertable, and protect the aquifer water resources.  However, the advancing speed of coalface is not the greater the better. The research shows that when the advancing speed is too high, the stress concentration in coal-rock mass is obvious, and a lot of capacity is accumulated at the same time (Wang et al. 2006;Xie et al. 2007;Xu et al. 2019;Li et al. 2021). Once the energy in coal-rock mass is released suddenly, it may lead to some dynamic disasters, such as rock burst, coal and gas outburst, mine earthquake, threatening the safety of miners (Zhang et al. 2021). Therefore, only an appropriate advancing speed of coalface can achieve the protection of shallow water resources under the premise of ensuring the safety of mine production.

Advancing speed optimization
Based on the mechanism analysis of advancing speed impact on overburden permeability and watertable of shallow aquifer, the advancing speed of panel 21103 in Yili No.4 Coal Mine, Xinjiang, China is optimized to improve the production efficiency of panel coalface, reduce the water inflow into coalface and protect the water resourcein shallow aquifer. Before the optimization of advancing speed, the allowable drawdown of shallow aquifer watertable in Yili No.4 Coal Mine should be determined. Liu et al. pointed out that the maximum watertable drawdown in Taklimakan Desert shall be no more than À3 m or else falling groundwater cannot prevent surface vegetative covers from drying up due to the lack of water (Liu et al. 2013). Taklimakan Desert is located in Xinjiang, China, which is near to Yili No.4 Coal Mine and has a certain similarity in surface vegetation. What's more, the desert area is arid climate with less precipitation and surface runoff. The À3 m decline threshold of watertable in Taklimakan Desert has some surplus for vegetation growth in Yili No.4 Coal Mine area. Here, this value is taken as the ultimate watertable drawdown which is able to  Figure 14. Also, the changes of watertable in gravel aquifer before and after the panel 21103 coalface advances through the borehole (from June 25, 2017 to July 19, 2017) are shown in Figure 15.
It can be seen from Figure 15 that, when coalface is at À60 m distance from the borehole (on June 25, 2017), aquifer watertable drawdown starts to appear. This phenomenon can be derived from the front disturbance induced by longwall coal mining and by reason that aquifer water flows into the subsidence basin (as shown in Figure  14), resulting in the decrease of the watertable. When the coalface arrives at the position just below the borehole (at 0 m distance from the borehole on July 3, 2017), aquifer watertable falls by À3.0 m; and when coalface advances 25 m beyond the borehole (on July 6, 2017), watertable further falls to À4.8 m and reaches the largest vertical drawdown. Although coal mining process continues after that, aquifer watertable gradually recovers on account of recharge from surface runoffs and melt water, and then returns to the peak value and maintains with coalface advancing 130 m beyond the borehole (on July 19, 2017). In the end, aquifer watertable is just À0.35 m lower than its original level. This shows that when the advancing speed of panel 21103 is 8 m/d, the shallow mudstone aquiclude has good water-resisting ability, which plays a good role in protecting the water resources of overlying gravel aquifer. It also means that the hydrological circumstance provided by shallow aquifer can support the regular growth of surface vegetation at the optimized advancing speed. After advancing speed of panel 21103 is optimized, artificial grass planting can be replaced by natural growth, which is demonstrated by the comparison between both advancing speeds shown in Figure 16.

Discussion
The effect of advancing speed on movement and permeability of overburden is actually an issue of time effect. Reflecting the different advancing speeds in numerical simulation is one of the problems that experts and scholars have been trying to solve. In order to reflect the difference of advancing speed in numerical simulation, the following two methods are mainly adopted by experts. Setting the same operation timesteps under the different excavation step spacing, or setting the different operation time-steps under the same excavation step spacing (Wang et al. 2006;Xie et al. 2007;Xu et al. 2019;Li et al. 2021). At the same time, combined with the verification of the numerical model, the above two methods can be used to achieve the similar simulation of advancing speed of coalface. Peer experts usually use one of the two methods to set the advancing speed of coalface when studying similar problems (Zhang et al. 2009;Tsesarsky 2012;Fan et al. 2018;Li et al., 2019;Guo et al. 2019;Wang et al. 2020). However, the above two methods can not achieve the complete agreement between the numerical model and the actual situation, which should be further improved and optimized in the future research.
According to this work, the break distance of overburden is shorter and the fracture development is more adequate after the advancing speed of coalface is reduced. However, the overburden breaking and fracture development have a time process. When the advancing speed of coalface decreases suddenly, the overlying rock does not produce short-distance broken rock blocks immediately. At the same time, the development of overburden fractures and the flow of shallow water from the aquifer to the coalface require a certain time. This makes that the increase of water inflow in coalface lags behind the decrease of advancing speed in Figure 2.

Conclusions
For weakly cemented strata, there is a butterfly-shaped permeability-increasing area that is comprised of a trapezoid and double V-shaped subareas because of longwall coal mining disturbance. The permeability-increasing subarea is symmetric about the centreline of the panel. The dominant seepage channels are located at the position corresponding to both ends of the panel coalface. With advancing speed slowing down, aquiclude would gradually lose its original integrity and experience 'effective water-resisting thickness getting thinfracture developmentfracture penetration' successively. Permeability curves of aquiclude are of double peaks. With the increase of panel advancing speed, the maximum permeability of aquifer exponentially decreases.
Under the circumstance of fixed head, the contour of subsurface watertable drawdown is of oval profile. With the increase of panel advancing speed, the maximum vertical drawdown shows a significant decrease obeying an exponential function y ¼ -5.68e À0.21x -1.76 (where y is the maximum vertical drawdown of aquifer watertable, m; x is advancing speeds, m/d; and the goodness of fit equals to 0.99). The fitting results show that, if taking À2 m vertical drawdown as an indicator, the lateral influence range linearly decreases, obeying y ¼ -8.75x þ 156.10 (where y is the lateral influence range, m; x is advancing speeds, m/d; and the goodness of fit equals to 0.98).
After underground coal mining, overburden strata collapses and forms stable voussoir beam structure. The voussoir beam is subjected to horizontal compressive force, which has a bearing capacity on the overlying strata. Strong supporting capacity of voussoir beam at a high advancing speed can control overlying strata to produce slight mining fractures, and keep aquiclude in a good integrity. Within a certain range, the higher the coalface advanced, the shorter time the cracks propagated, the faster the cracks closed, and the more gently the surface subsided. In addition, for weakly cemented strata, the rich clay minerals such as montmorillonite and kaolinite can be expanded and argillated with water, which can block some tiny mining cracks and has a good closing effect on the slight cracks produced at a high advancing speed. This is the reason why the higher the advancing speed, the lower the overall permeability of overburden, and the less the aquifer watertable drops to some extent. Therefore, within a certain range, improving the advancing speed of panel coalface can reduce the development of overlying strata fractures and protect the aquifer water resources.
By taking the ultimate watertable drawdown that is still able to sustain the regular growth of surface vegetation as a constraint condition, the optimal advancing speed of Panel 21103 is determined as 8 m/d. In such a scenario, aquifer watertable would witness a temporary fall, but then recover due to recharge originating from surface runoffs and melt water. The optimized advancing speed contributes to the shift from artificial planting to natural growth of the surface vegetation in mining area and as a result, local ecological environment can be protected.