Numerical inversion of stress evolution of mining strata and stress path equivalent conversion method
By simulating the coal seam excavation process in continuous media analysis software, and using the Double Yield material model and equivalent conversion method, the problem of simulating the collapse, accumulation and compaction of the roof in the goaf in continuous media software was solved. This enabled accurate simulation and conversion of stress paths throughout the entire coal seam mining process, providing a foundation for coal and rock mechanical property tests under mining-induced stress paths.
Patent Information
- Application Number
- CN202210553530.3
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2022-05-20
- Publication Date
- 2026-02-06
- Estimated Expiration
- 2042-05-20
AI Technical Summary
Existing continuous medium software is unable to effectively simulate the collapse, accumulation, and compaction processes of the goaf roof, resulting in inaccurate stress path simulation during coal seam mining and an inability to accurately reflect the actual mining stress environment.
By performing secondary development in continuous media analysis software, a numerical model was established using FLAC3D, monitoring points were set, and roof collapse and gangue accumulation during coal seam excavation were simulated. The Double Yield material model was used to simulate gangue compaction bearing capacity, and the mining-induced stress path was converted into the laboratory loading and unloading stress path through the equivalent transformation method.
It achieves accurate simulation of roof and floor stress and equivalent transformation of stress path throughout the entire coal seam mining process, provides a basis for testing coal and rock mechanical properties under mining-induced stress path, and solves the problem of continuous medium software in simulating stress evolution in goaf.
Smart Images
Figure CN114925524B_ABST
Abstract
Description
TECHNICAL FIELD
[0001] The present application relates to a kind of mining rock stress evolution numerical inversion and stress path equivalent conversion method, especially suitable for coal mining technical field. BACKGROUND
[0002] The mechanical properties such as damage, failure and seepage of coal rock are closely related to the stress path it undergoes, and the physical and mechanical properties of coal rock under different stress paths are necessarily different. The stress path in traditional tests is mainly for static environment, while coal mining is a dynamic stress environment. Therefore, it is necessary to study the stress evolution law of surrounding rock in the whole process of coal mining, and then convert it into laboratory stress path to carry out coal rock mechanical property test under mining stress path, and obtain the characteristics of mining coal rock damage, failure, etc. to provide support for guiding field mining practice.
[0003] At present, the research on the stress evolution law of the surrounding rock in the advanced section during coal mining has been relatively mature, while the research on the stress evolution law in the goaf after coal mining is still in the exploratory stage. Due to the concealment and complexity of the goaf, it is difficult to study it by field measurement method. Therefore, it is a feasible method to study the stress recovery characteristics of the goaf by using numerical simulation software.
[0004] The present application uses continuous medium analysis software widely used in engineering field, simulates the stress evolution law of roof and floor in the whole process of coal mining by secondary development of continuous medium analysis software, and converts it into laboratory loading and unloading stress path to provide basis for studying the mechanical properties of coal rock under mining stress path. Compared with the existing research, the main innovations of the present application are as follows:
[0005] (1) Roof caving process simulation
[0006] The blocks in continuous medium software cannot be separated from each other, so it is impossible to directly simulate the caving and accumulation process of goaf roof. In the existing research, the roof caving is simulated by directly deleting the unit bodies in the corresponding height range, for example, in the document "Mining Response Inversion Based on Compaction Theory of Goaf", an average caving height is calculated by an empirical formula, and then the unit bodies in the corresponding height range are directly deleted to simulate the roof caving. However, in the actual coal mining process, the roof at different positions of the working face has different lithology, support conditions and stress states, so the caving height is also different, that is, the average caving height does not reflect the characteristics of roof caving with mining in the actual mining process.
[0007] (2) Simulation of caving gangue accumulation process
[0008] The simulation of the accumulation process of the caved gangue in the goaf is mostly through filling the unit bodies of different material models, for example, in the document "Research on the coupling analysis method of mining stress and goaf compaction bearing", the accumulation process of the gangue is simulated by filling the unit bodies of the Double Yield material model in the goaf. However, in fact, not all caved gangues can form support to the roof at different positions, and the accumulation height of the caved gangue at some positions is less than the caving height, in which case the caved gangue does not contact the roof and cannot form effective support to the roof. According to the existing research results, the roof at the boundary of the goaf is not fully caved, and the gangue cannot form support to the roof; and the roof in the middle region of the goaf is fully caved, and the gangue is easy to form effective support to the roof. Therefore, if only the method of filling all the space in the goaf is used, the actual accumulation state of the gangue in the goaf cannot be reflected.
[0009] (3) Simulation of the compaction bearing process of the accumulated gangue
[0010] For the simulation of the mechanical properties of the compaction bearing process of the accumulated gangue in the goaf, the existing research shows that the mechanical properties can be represented by the Salamon relationship, and the Double Yield material model is generally used to simulate the compaction bearing properties of the accumulated gangue. The model can better simulate the plastic hardening characteristics of the material, which is similar to the compaction process of the accumulated gangue in the goaf. However, in most studies, the material parameters of the Double Yield model are calibrated by the trial-and-error method and the iterative method, for example, in the document "Numerical simulation method of mine pressure in the whole process of longwall mining". Thus, not only is the work large, but also the simulation effect varies greatly under different geological conditions, and it is not universally applicable. SUMMARY
[0011] Technical purpose: In view of the shortcomings of the prior art, a mining rock stress evolution numerical inversion and stress path equivalent conversion method is provided, which can overcome the difficulty of the continuous medium software in simulating the caving, accumulation, compaction and bearing processes of the roof in the goaf, realize the numerical simulation of the whole coal mining process, obtain the stress evolution law of the roof and floor in the whole process, and equivalent convert the test loading and unloading stress path, thereby providing a mining rock stress evolution numerical inversion and stress path equivalent conversion method for studying the coal and rock mechanical properties under the mining stress path.
[0012] Technical scheme: To achieve the above object, the mining rock stress evolution numerical inversion and stress path equivalent conversion method of the application firstly, according to the actual mining geological conditions, the corresponding numerical model is established by using computer software, and the monitoring points are set in the different layers of the roof and floor of the numerical model; Then, the coal seam in the model is simulated to dig, and the roof is collapsed and the gangue is accumulated in real time with the digging, so that the whole process of the coal seam, the roof and floor changing with the mining is inverted, and the stress change law of the roof and floor in the whole process of coal mining is obtained, wherein the vertical stress and horizontal stress of the roof present opposite change law as a whole, and the vertical stress and horizontal stress of the floor present the same change law; Secondly, according to the synchronous change law of vertical and horizontal stress, the stress path characteristics of mining roof and floor are summarized, wherein the stress path of the roof presents 'positive S shape', and the stress path of the floor presents'reverse S shape'; Finally, through the equivalent conversion method, the vertical stress and horizontal stress are equivalent to the axial pressure and confining pressure, and the stress evolution process of the roof and floor is equivalent to the triaxial loading and unloading stress path in the test, so as to realize the equivalent conversion of the mining stress path.
[0013] The specific steps are:
[0014] Step one, the coal rock sample of the stratum in the region to be simulated is obtained by field drilling, the mechanical property test is carried out, the basic mechanical parameters of the coal seam and the roof and floor rock stratum are obtained, according to the actual geological conditions, the geological model is established by using FLAC3D, a plurality of monitoring points are set in the different layers of the roof and floor of the geological model, and the stratum of the coal seam in the geological model is composed of a plurality of basic units;
[0015] Step two, the simulation of the coal seam in the geological model is excavated, the simulation excavation strategy is set according to the need, that is, the roof collapse program is executed after every a meter of excavation is calculated b time steps; That is, through the excavation strategy, the volume strain ε v and principal stress of each basic unit are obtained, the damage value of each basic unit is calculated; When the damage value is greater than the preset critical damage value D L , it is considered that the basic unit has been damaged, the material model of the basic unit is replaced by null, and it is defined as the roof collapse zone Cave group; In this way, all the basic units in the geological model are judged;
[0016] Step three, the top judgment of the goaf is carried out: the z axis coordinate Z max of the maximum height of the collapse zone at the basic unit coordinates (x, y) of all the collapse Cave group is calculated, the z axis coordinate Z min of the coal seam floor is subtracted, the collapse zone rock stratum thickness ∑h i is obtained, the collapse zone rock stratum thickness ∑h i is multiplied by the dilatancy coefficient B of the rock stratum to obtain the accumulation height h of the collapsed rockc ; compare the height h of the caving zone k and the height h of the pile c , when h c ≥ h k , it means that the rock of the caving zone roof has been able to form a compact support for the roof at this time, at this time the supported roof caving zone is replaced by the Double Yield model, which is defined as Fillgroup; if h c < h k , it means that the caving gangue has not formed a support for the roof, and no operation is performed; the caving, roof joining and compaction state of the goaf roof is repeatedly judged every time the b step is calculated;
[0017] Step four, carry out gangue compaction bearing judgment: first, obtain the basic mechanical parameters such as the density p, tensile strength s t , bulk modulus K of the accumulated gangue collected in the simulation area through laboratory coal rock sample mechanical property test; then determine the coefficient R, cap pressure p c , plastic strain increment e p according to the actual mining geological conditions and the mechanical properties of the broken gangue in the goaf Salamon constitutive relationship; assign the above parameters to the basic unit defined as Fill group through fish language;
[0018] Step five, continue to perform mechanical calculation on the model that has not reached mechanical equilibrium, and judge whether the maximum unbalanced force of the geological model is less than 1 -6 ; if the geological model reaches equilibrium, judge whether the geological model is excavated; if not, continue to execute the next excavation stage, and repeatedly perform roof caving calculation, gangue accumulation calculation and compaction bearing calculation;
[0019] Step six, complete the excavation of the geological model;
[0020] Step seven, collect the stress change data of the monitoring points arranged in the roof and floor of different layers in the whole process of coal seam mining, and draw the evolution law of the vertical stress s z and the horizontal stress s x according to the collected stress change data, so as to deduce that the stress path of the mining roof presents a positive "S" shape, and the stress path of the mining floor presents a negative "S" shape;
[0021] Step eight, equivalent the axial pressure in the triaxial test of the strata coal rock sample to the vertical stress s z , and equivalent the confining pressure to the horizontal stress s x , convert the stress path of the mining roof and floor into the loading and unloading stress path in the laboratory, and carry out the mechanical property test of the mining coal rock.
[0022] Further, the specific steps of executing the roof caving procedure are as follows:
[0023] Through the continuous medium numerical software, the fish language built-in is used to compile and realize:
[0024] a1, through the fish language compilation algorithm, all basic units in the geological model are traversed, and the number id, coordinates (x, y, z), principal stress (σ1, σ2, σ3), volume strain ε of each basic unit are obtained v , the axial pressure σ1 is equivalent to the vertical stress, and the confining pressure σ2 and σ3 are equivalent to the horizontal stress;
[0025] a2, using the obtained principal stress (σ1, σ2, σ3), volume strain ε v , the damage degree D of each basic unit is calculated;
[0026] a3, when the damage degree D of the basic unit is greater than or equal to the preset critical caving threshold D L , it is considered that the basic unit has been damaged, and the material model of the basic unit is replaced by a null model, that is, the roof rock layer is simulated to caving, and the unit body is defined as Cave group; if the damage degree D is less than the critical caving threshold D L , it is considered that the unit body has not been damaged, and no treatment is performed;
[0027] a4, the above three steps are cycled for all basic units until the judgment of all basic units in the geological model is completed.
[0028] Further, the damage degree D in step a2 is calculated by the following formula:
[0029] ;
[0030] The critical caving threshold D in step a3 is calculated by the following formula: L :
[0031]
[0032] In the formula: ε V e , ε V p respectively represent the elastic body strain and the plastic body strain, E represents the elastic modulus, v represents the Poisson's ratio, σ r , σ f respectively represent the residual strength and the peak strength.
[0033] Further, the top of the gob is judged by calculating the gangue accumulation process, and the algorithm steps compiled by the fish language built-in in the continuous medium numerical software are as follows:
[0034] b1. Using the FILE programming language, implement an algorithm to obtain the maximum z-coordinate Z at each coordinate point in the x and y planes of all basic units in the Cave group. max The minimum coordinate Z at that coordinate location min That is, the bottom plate of the coal seam;
[0035] b2. Obtain the coordinates Z max Z min Substitute into the formula: ∑h i =Z max -Z min -M, the total thickness ∑h of the collapsed rock layer can be calculated. i Let the coefficient of fragmentation be B, and then calculate the accumulation height h at that coordinate. c =B∑h i Collapse height h k = Z max -Z min ;
[0036] b3. When the roof accumulation height of the goaf is h c Greater than or equal to the collapse height h k When the accumulated gangue reaches a certain height, it indicates that the accumulated gangue has come into contact with the roof and formed an effective support for the roof. In FLAC3D software, the material model along the entire z-direction at this coordinate within the goaf is changed to a Double Yield model, and the basic element is defined as a Fill group. c Less than the collapse height h k At this time, the accumulated gangue has not yet formed an effective support for the roof and does not provide any force. Therefore, it has no impact on the overall model calculation effect and no treatment is required for the goaf in this case.
[0037] b4. Repeat steps b1 to b3 above to iterate and judge the coordinates in the entire goaf space, thereby realizing the numerical simulation of the gangue accumulation process during the simulation of mining in the entire geological model.
[0038] Furthermore, the numerical implementation method for the compaction bearing capacity of stockpiled gangue is as follows:
[0039] c1. Plastic hardening parameters for the Double Yield material model simulating accumulated gangue: coefficient R, cap pressure p c Plastic strain increment ε p Perform calibration;
[0040] c2. Using the fish language, mechanical parameters are assigned to the basic elements in the Fill group to simulate the compaction bearing mechanical characteristics of gangue in the goaf.
[0041] Furthermore, step c1 specifically includes:
[0042] c11, sampling the goaf accumulated broken gangue on site in the simulation area, testing the relevant mechanical properties of the accumulated broken gangue, and obtaining the basic mechanical parameters of the accumulated broken gangue, including density p, cohesive force c, internal friction angle f, tensile strength s t, bulk modulus K, and shear modulus G;
[0043] c12, using a uniaxial compression experimental device of the accumulated broken gangue, obtaining stress and strain data in the compression process of the accumulated broken gangue;
[0044] c13, the coefficient R in the Double Yield material model controls the ratio of the axial pressure to the confining pressure;
[0045] c14, using the Salamon constitutive relation to represent the compression bearing characteristics of the accumulated broken gangue, substituting the obtained stress and strain data into the calculation formula of the cap pressure p c and the plastic strain increment e p to obtain the corresponding values.
[0046] The coefficient R in step c13 is calculated by the following formula:
[0047] ;
[0048] The cap pressure in step c14 is calculated by the following formula:
[0049] ;
[0050] The plastic strain increment in step c14 is calculated by the following formula:
[0051]
[0052] wherein s represents the axial pressure, s
[0053] Further, the mining stress path equivalent conversion method refers to using the pressures in three directions in the triaxial test instrument to equivalently replace the stresses monitored by the element in the numerical model, i.e. s z , s y , and s x , equivalently converting the stress evolution law obtained in the process of coal mining into a mining stress path, and restoring the actual mining stress environment through the triaxial test instrument.
[0054] Beneficial effects:
[0055] 1) The method can solve the problem that the roof caving, accumulation, compaction and bearing process in the coal seam mining process cannot be simulated in the continuous medium software.
[0056] 2) The roof and floor mining stress path in the whole process of coal seam mining can be obtained, and the foundation for developing the coal rock mechanics test of mining coal is provided.
[0057] By defining the caving damage threshold of the roof, then judging the damage degree of the roof through the loop traversal algorithm, and comparing with the caving threshold, the simulation of the roof caving process in the mining process is realized, the corresponding relationship between the Double Yield material model related parameters and the Salamon classic model is established, and the accurate calibration of the Double Yield model material parameters is realized. Thus, the numerical model of the inversion of the roof caving, accumulation, compaction and bearing process of the goaf is developed, and the stress evolution law of the roof and floor in the whole process of coal seam mining is fully reflected. The whole process of mining stress evolution is equivalent to the laboratory loading and unloading stress path, and the actual mining mechanics environment is simulated. BRIEF DESCRIPTION OF DRAWINGS
[0058] Figure 1 The figure is the flowchart of the mining rock stress evolution numerical inversion and stress path equivalent conversion method of the application.
[0059] Figure 2 The figure is the schematic diagram of the monitoring point arrangement position in different mining stages.
[0060] Figure 3 The figure is the roof vertical / horizontal stress evolution law.
[0061] Figure 4 The figure is the schematic diagram of the roof caving zone vertical / horizontal stress synchronous change.
[0062] Figure 5 The figure is the schematic diagram of the roof fracture zone vertical / horizontal stress synchronous change.
[0063] Figure 6 The figure is the schematic diagram of the roof bending subsidence zone vertical / horizontal stress synchronous change.
[0064] Figure 7 The figure is the schematic diagram of the mining roof stress path.
[0065] Figure 8 The figure is the schematic diagram of the floor vertical / horizontal stress evolution law.
[0066] Figure 9 The figure is the schematic diagram of the floor mining damage zone vertical / horizontal stress synchronous change.
[0067] Figure 10This is a schematic diagram showing the synchronous change of vertical and horizontal stresses in the bending deformation zone of the base plate.
[0068] Figure 11 This is a schematic diagram of the stress path in the mining base plate.
[0069] In the figure: 1-Monitoring point of roof collapse zone, 2-Monitoring point of roof fracture zone, 3-Monitoring point of roof bending and subsidence zone, 4-Monitoring point of floor mining-induced damage zone, 5-Monitoring point of floor bending and deformation zone, 6-Roof, 7-Coal seam, 8-Floor, 9-Goaf, 10-Roof bending and subsidence zone, 11-Roof fracture zone, 12-Roof collapse zone, 13-Floor mining-induced damage zone, 14-Floor bending and deformation zone. Detailed Implementation
[0070] The embodiments of the present invention will be further described below with reference to the accompanying drawings:
[0071] like Figure 1 As shown, the present invention provides a method for numerical inversion of stress evolution in mining-induced rock strata and equivalent transformation of stress paths, comprising the following steps:
[0072] S1. Drilling in the field to obtain coal and rock samples from the strata, and conducting mechanical property tests to obtain the basic mechanical parameters of the coal seam and the roof and floor strata;
[0073] S2. Based on the actual geological conditions of the mining operation, a numerical model is established using continuous media software. Monitoring points are arranged at different levels in the roof and floor of the coal seam, including monitoring points 1 (roof collapse zone), 2 (roof fracture zone), and 3 (roof bending and subsidence zone) at different locations in the roof 6, and monitoring points 4 (floor mining-induced failure zone) and 5 (floor bending and deformation zone) at different locations in the floor 8 below the coal seam 7. Figure 2 As shown;
[0074] S3. Simulate coal seam mining in the numerical model, and execute the roof collapse calculation program, gangue accumulation calculation program, and compaction bearing capacity calculation program written in the fish language in a loop, and collect stress data of all monitoring points in real time throughout the process.
[0075] S4. After the coal seam is mined out, the stress evolution paths of the roof and floor are converted into laboratory loading and unloading stress paths using the equivalent conversion method. Information on the roof bending and subsidence zone 10, roof fracture zone 11, roof collapse zone 12, floor mining-induced failure zone 13, and floor bending deformation zone 14 is obtained through analysis, thus acquiring various information on stress changes. Figure 2 In this context, A represents the model state of the original rock stress. Figure 2 In this model, B represents the stress-increased state. Figure 2 In the model, C represents the stress reduction state. Figure 2 In this context, D represents the model state of stress recovery.
[0076] Firstly, according to the actual mining geological conditions, the corresponding numerical model is established by using computer software, and monitoring points are set at different layers of the roof and floor of the numerical model; then the coal seam in the model is simulated to be mined, and the roof collapses in real time and the gangue is accumulated with the mining, thereby the whole process of the change of the coal seam, roof and floor with the mining is inversed, and the stress change law of the roof and floor in the whole process of the coal seam mining is obtained, wherein the vertical stress and the horizontal stress of the roof present opposite change law as a whole, and the vertical stress and the horizontal stress of the floor present the same change law; secondly, the stress path characteristics of the mining roof and floor are summarized according to the synchronous change law of the vertical and horizontal stresses, as shown in Figure 4 、 Figure 5 and Figure 6 , wherein the stress path of the roof presents a "positive S shape", and the stress path of the floor presents a "reverse S shape", as shown in Figure 7 ; finally, the vertical stress and the horizontal stress are equivalent to the axial pressure and the confining pressure by the equivalent conversion method, and the stress evolution process of the roof and floor is equivalent to the triaxial loading and unloading stress path in the test, so as to realize the equivalent conversion of the mining stress path, as shown in Figure 3 .
[0077] The evolution law of the vertical / horizontal stress of the floor is shown in Figure 8 , the synchronous change of the vertical / horizontal stress of the mining damage zone of the floor is shown in Figure 9 , the synchronous change of the vertical / horizontal stress of the bending deformation zone of the floor is shown in Figure 10 , and the stress path of the mining floor is shown in Figure 11 . The above changes all represent different deformation and damage zones, and the damage form of the rock in each subzone is different, and the stress evolution may also be different. Therefore, the corresponding monitoring points are selected in each subzone by the method, the stress change is analyzed, and the stress evolution law of the point is used to represent the stress change of the roof and floor at different layers.
[0078] As shown in Figure 1 , the specific steps include the numerical implementation method of the roof collapse process, the numerical implementation method of the gangue compaction process in the goaf, the numerical implementation method of the compaction and bearing process of the accumulated gangue, and the equivalent conversion method of the mining stress path:
[0079] 1) Numerical implementation method of roof collapse
[0080] Although the block in the continuous medium software cannot collapse by itself, it can be replaced by an empty model to simulate the roof collapse when the block meets the collapse condition by setting the judgment condition. The damage of rock is often associated with damage, so the damage value of rock is used to measure the damage of rock, and then the collapse of the roof is judged. According to the plastic theory, the general definition of damage is: D=1-A0 / A, wherein: A is the bearing area of the complete (no damage) rock; A0 is the effective bearing area of the rock after damage.
[0081] The deformation of rock after damage in the process of stress can be divided into elastic and plastic deformation. Among them, plastic deformation will lead to plastic body strain ε v p Increase, so A= A0(1+ε v p ), the damage variable can be expressed as: D=ε v p / (1+ε v p ). Since ε v p <<1, in order to simplify the calculation, the damage variable is approximately equal to the plastic body strain, and the plastic body strain is equal to the total volume strain minus the elastic body strain, so:
[0082] (1)
[0083] Where: σ1, σ2, σ3 are the maximum, intermediate and minimum principal stress, v is the Poisson's ratio, ε v p is the plastic body strain of rock, ε v e is the elastic body strain of rock, and E is the elastic modulus.
[0084] Therefore, by compiling algorithms in fish language to traverse all elements in the model, the volume strain ε v and principal stress of each element are obtained, and the damage degree D of each element is calculated. When D is greater than or equal to the critical collapse threshold D L , it is considered that the element has been damaged, and the material model of the element is replaced by an empty model, that is, the roof rock layer is simulated to collapse, and the element is defined as Cave group. If the damage degree D is less than the critical collapse threshold D L , it is considered that the element has not been damaged, and no treatment is performed.
[0085] 2) Numerical implementation method of gangue accumulation
[0086] Not all collapsed gangue in the goaf can form support to the roof, and the bulking coefficient B of the collapsed gangue needs to be considered, and the accumulation height h c=∑h i If the calculated accumulation height h c is greater than the caving height h k =∑h i +M, the accumulated gangue forms an effective support to the roof, and the simulation of the gangue accumulation is realized by changing the material model of the goaf space; if h c <h k , it indicates that the accumulated gangue has not contacted the roof, and this case is not processed.
[0087] According to the method, the maximum z coordinate Z max at each coordinate in the x, y plane in the Cave group is obtained by the fish language algorithm. min The minimum coordinate (coal seam floor) corresponding to the coordinate is Z i , and the thickness of the caving strata is ∑h max =Z min -Z c -M, assuming that the crushing expansion coefficient has been given, the accumulation height h k and the caving height h m at the coordinate can be calculated; then it is judged whether it meets the support condition by the above method, and if it meets the support condition, the material model of the entire z direction within the goaf at the coordinate is changed to the Double Yield model. The coordinates within the entire goaf space are judged by the loop algorithm compiled by the fish language, thereby realizing the numerical simulation of the gangue accumulation process.
[0088] 3) Numerical implementation method of the bearing process of the accumulated gangue
[0089] The process of the goaf gangue compaction bearing is similar to the axial compaction experiment of the broken rock in the laboratory, and therefore, the compaction mechanical properties of the broken rock obtained in the laboratory can be used to characterize the mechanical behavior of the goaf gangue. The caving and accumulated gangue in the goaf gradually compacts in the process of the bending subsidence of the overlying roof, as shown in Figure 4 , and the stress also gradually increases with the increase of the vertical strain, but the change process is not a simple linear increase process, but a nonlinear cumulative process. Existing researches show that the Salamon relationship σ=E0 / (1-ε / ε m ) can better characterize the compaction bearing mechanical properties of the accumulated gangue. Therefore, how to establish the relationship between the material parameters of the accumulated gangue in the numerical model and the Salamon parameters is the key to realize the simulation of the compaction bearing process of the gangue.
[0090] Take the Double Yield model as an example, which can be used to characterize the irreversible compression behavior of materials after being stressed. By analyzing the correspondence between the relevant parameters in the Double Yield material model and the parameters in the Salamon classic model of the compaction process of broken rocks, it is accurately calibrated. The material parameters that need to be determined in the Double Yield material model mainly include density ρ, cohesive force c, internal friction angle φ, tensile strength σt, bulk modulus K, shear modulus G, coefficient R, plastic strain increment ε p , cap pressure p c . Among them, the density ρ is a basic parameter of the material, which can be selected according to the actual lithology. The cohesive force c, the internal friction angle φ, and the tensile strength σ t can be obtained through basic mechanical property tests of coal rock samples. The bulk modulus K and the shear modulus G can be determined through their relationship with the elastic modulus E and the Poisson's ratio v, K=E / 3(1-2v), G=E / 2(1+v).
[0091] The key is how to determine the coefficient R, the cap pressure p c , and the plastic strain increment ε p . The coefficient R is the ratio of the elastic bulk modulus to the plastic bulk modulus of the material, which mainly determines the slope of the stress-strain unloading curve when the material enters the plastic stage. However, in actual production, the boundary of the goaf is often not filled, so the pressure on the surrounding accumulated gangue is limited, and σ3 can be considered as 0. According to the relationship between the coefficient R and the axial pressure and the confining pressure:
[0092] (2)
[0093] The coefficient R can be determined as R=1.5.
[0094] In the Double Yield model, the introduction of the volume yield surface (cap) realizes the irreversible volumetric deformation after being stressed, where the volume yield surface is defined by p c ("cap pressure"), which satisfies the following relationship with the material properties and the current stress state:
[0095] (3)
[0096] And the plastic strain increment ε p satisfies the following relationship with the current vertical strain of the material:
[0097] (4)
[0098] In this way, the cap pressure p c and the plastic strain increment ε p can be determined through equations (3) and (4).
[0099] According to this method, if the actual geological conditions of a coal seam are known, the Salamon relationship satisfied during the compression of its fractured roof can be obtained through relevant experiments, thus yielding the initial elastic modulus E0 and the maximum compressive strain ε. m Substituting these equations into (3) and (4) will determine p. c ε p Parameters. By assigning corresponding material parameters to the stockpiled gangue material model during the simulation process, the compaction and bearing process of the specified stockpiled gangue can be accurately simulated.
[0100] 4) Equivalent transformation method for dynamic stress path
[0101] The above three methods can be used to numerically simulate the entire process of roof collapse, accumulation, compaction, and load-bearing during coal seam mining. Combined with actual geological conditions of coal seam mining, corresponding numerical models can be established to simulate the entire excavation process, and stress change data at any location on the roof and floor of the coal seam can be obtained through monitoring.
[0102] Taking the coal seam roof as an example, a monitoring point is set 20 m above the roof in the middle of the working face. The vertical stress σ is then statistically analyzed throughout the entire process from -150 m to 150 m from this point (from the advance section to the post-mining stage). z Horizontal stress σ x σ y Data. This allows us to obtain the stress path of the roof during mining. In a true triaxial test in the laboratory, the stress of the element in the equivalent numerical model can be obtained using pressure in three directions, i.e., σ1 = σ z σ²=σ y σ3=σ x If the laboratory conditions are a pseudo-triaxial test (σ2=σ3), the vertical stress σ can be analyzed. z With horizontal stress σ x The variation law, with axial compression σ1 equivalent to σ z The confining pressure σ3 is equivalent to σ x .
[0103] Similarly, the evolution law of mining-induced stress in the floor can be obtained separately throughout the entire coal seam mining process, such as... Figure 3 As shown, this method converts the evolution law of mining-induced stress into a laboratory mining-induced stress path, simulating the actual mining disturbance stress environment.
[0104] Example 1:
[0105] The first step is to establish a reasonable calculation model based on the actual geological conditions and set up monitoring points at different layers of the top and bottom slabs.
[0106] The second step involves calculating b steps (e.g., 500 steps) and executing the roof collapse procedure every am (e.g., 5 m) of excavation. This means using an algorithm to traverse all elements in the model and obtain the volumetric strain ε of each element. v Based on the principal stresses, the damage value of each element is calculated. When it exceeds the critical damage value D... L If the element is found to be damaged, its material model is replaced with null, and it is defined as a Cave group. This process is repeated for all elements in the model.
[0107] The third step is to execute the goaf roof contact judgment procedure. The algorithm obtains the z-axis coordinate Z of the maximum height of the collapse zone at coordinates (x, y) in the Cave group. max Subtract the z-axis coordinate of the coal seam floor. min The thickness of the strata in the collapse zone, ∑h, was obtained. i Multiplying by the rock strata's fragmentation coefficient B yields the accumulation height h of the collapsed rocks. c Compare the height h of the landslide zone. k and stacking height h c When h c ≥h k When this indicates that the collapsed rock has already provided support for the roof, the material model of all unit cells in the z-direction at (x, y) within the goaf space is replaced with a Double Yield model and defined as a Fill group. If h c <h k This indicates that the collapsed gangue has not yet provided support for the roof, and no further action should be taken.
[0108] The fourth step is to execute the gangue compaction and bearing capacity procedure. First, the mechanical properties of the accumulated gangue (ρ) and tensile strength (σ) are obtained through laboratory testing of coal and rock samples. t Basic mechanical parameters such as bulk modulus K; then, based on actual mining geological conditions and the Salamon constitutive relation, the coefficients R and cap pressure p are determined. c Plastic strain increment ε p The above parameters are assigned to the unit cells in the Fill group using the Fish language.
[0109] Fifth, continue with mechanical calculations to determine if the maximum unbalanced force of the model is less than 1. -6 If the model reaches equilibrium, determine whether the excavation is complete. If not, continue to the next excavation stage. The roof collapse calculation program, the gangue accumulation calculation program, and the compaction bearing capacity calculation program are executed cyclically.
[0110] Step 6: Model excavation completed.
[0111] The seventh step is to summarize the vertical stress σ z and the horizontal stress σ x of different layers of the roof and floor monitoring points in the whole process of coal mining, and find that the vertical stress and the horizontal stress of the roof present opposite change law, while the vertical stress and the horizontal stress of the floor present the same change law.
[0112] The eighth step is to draw the evolution law of the vertical stress σ z and the horizontal stress σ x , and find that the stress path of the mining roof presents a positive "S" shape, and the stress path of the mining floor presents a reverse "S" shape, as shown in Figure 6 .
[0113] The ninth step is to convert the stress path of the mining roof and floor into the stress path of the laboratory loading and unloading by using the axial pressure σ1 equivalent to the vertical stress σ z and the confining pressure σ3 equivalent to the horizontal stress σ x , and carry out the coal rock mechanical property test under mining.
Claims
1. A method for numerical inversion of stress evolution in mining-induced rock strata and equivalent transformation of stress paths, characterized in that: The specific steps are as follows: Step 1: Obtain strata coal and rock samples by field drilling in the area to be simulated, conduct mechanical property tests, obtain basic mechanical parameters of coal seam and roof and floor strata, and establish a geological model using FLAC3D based on actual geological conditions. Set multiple monitoring points at different layers of the roof and floor of the geological model. The strata of the coal seam in the geological model are composed of multiple basic units. Step 2: Simulate the excavation of the coal seam in the geological model. The simulation excavation strategy is set as needed to calculate step b after every 'a' meters of excavation, and then execute the roof collapse procedure. That is, the excavation strategy traverses all the basic units of each stratum in the geological model to obtain the volumetric strain ε of each basic unit. v Based on the principal stresses, the damage value of each basic element is calculated; when the damage value exceeds the preset critical damage value D... L When the basic unit is considered to have been destroyed, its material model is replaced with an empty model (null), and it is defined as the roof collapse zone (Cave group). This process is repeated for all basic units in the geological model. The specific steps for implementing the roof collapse procedure are as follows: This was implemented using the built-in fish language in continuous medium numerical software: a1. Using the FILE programming language, implement an algorithm to traverse all basic elements in the geological model and obtain the element ID, coordinates (x, y, z), principal stresses (σ1, σ2, σ3), and volumetric strain ε for each basic element. v σ1 is equivalent to vertical stress, and σ2 and σ3 are equivalent to horizontal stress; a2. Using the obtained principal stresses (σ1, σ2, σ3) and volumetric strain ε v Calculate the damage level D of each basic unit; ； The critical collapse threshold D in step a3 is calculated using the following formula. L : Where: ε V e ε V p Let E and σ represent the elastic strain and plastic strain, respectively, where E represents the elastic modulus, v represents Poisson's ratio, and σ represents the elastic strain and σ represents the plastic strain. r σ f Residual strength and peak strength, respectively; a3. When the damage level D of the basic unit is greater than or equal to the preset critical collapse threshold D L If the basic unit is considered to have failed, its material model is replaced with a null model, simulating the collapse of the roof strata, and the unit is defined as a Cave group; if the damage level D is less than the critical collapse threshold D... L At that time, it was assumed that the unit had not yet been damaged and no action was taken; a4. Repeat steps a1 to a3 above for all basic units until all basic units in the geological model have been judged. Step 3: Determine the roof contact in the goaf: Calculate the z-axis coordinate (Z) of the maximum height of the goaf zone at the basic unit coordinates (x, y) of all goaf groups. max Subtract the z-axis coordinate of the coal seam floor. min The thickness of the strata in the collapse zone, ∑h, was obtained. i Thickness of the caving zone rock strata ∑h i Multiplying by the rock strata's fragmentation coefficient B yields the accumulation height h of the collapsed rocks. c Compare the height h of the landslide zone k and stacking height h c When h c ≥h k When h indicates that the collapsed rock in the goaf caving zone has already formed a compacted support for the roof, the roof caving zone that has achieved support is replaced by a double yield model, which is defined as the Fillgroup; if h c <h k This indicates that the collapsed gangue has not yet provided support for the roof, so no further action is taken; after each cycle of calculation step b, the collapse, connection, and compaction status of the roof in the goaf are repeatedly assessed. Step 4: Determine the compaction bearing capacity of the gangue: First, obtain the density ρ and tensile strength σ of the accumulated gangue collected in the simulated area by testing the mechanical properties of coal and rock samples in the laboratory. t The bulk modulus K is a fundamental mechanical parameter; then, based on the actual mining geological conditions and the Salamon constitutive relation of the crushed gangue in the goaf, the coefficients R and cap pressure p are determined. c Plastic strain increment ε p The above parameters are assigned to the basic unit defined as a Fill group using the Fish language; Step 5: Continue mechanical calculations on the model that has not yet reached mechanical equilibrium to determine whether the maximum unbalanced force of the geological model is less than 1. -6 ; If the geological model reaches equilibrium, determine whether the geological model has been excavated. If the excavation has not been completed, continue to the next excavation stage and cycle through the calculations of roof collapse, gangue accumulation, and compaction bearing capacity. Step Six: Complete the geological model excavation; Step 7: Collect stress change data from monitoring points located in different layers of the roof and floor throughout the entire coal seam mining process, and plot the vertical stress σ based on the collected stress change data. z With horizontal stress σ x The evolution law is used to deduce that the stress path of the mining roof is generally positive "S" shaped, and the stress path of the mining floor is generally negative "S" shaped. Step 8: Equivalent the axial compression of the formation coal and rock samples in the triaxial test to the vertical stress σ. z The confining pressure is equivalent to the horizontal stress σ. x The stress paths of the mining-induced roof and floor plates are converted into laboratory loading and unloading stress paths to conduct tests on the mechanical properties of the mining-induced coal and rock.
2. The numerical inversion method for stress evolution of mining-induced strata and equivalent transformation of stress paths as described in claim 1, characterized in that... By calculating the rockfill process, the roof contact of the goaf can be determined. The algorithm steps, developed using the built-in FISH language in the continuous medium numerical software, are as follows: b1. Using the FILE programming language, implement an algorithm to obtain the z-axis coordinate of the maximum height at each coordinate in the x and y planes of all basic units in the Cave group. max The z-axis coordinate of the coal seam floor at that coordinate is Z. min That is, the bottom plate of the coal seam; b2. Obtain the z-axis coordinate of the maximum height. max Z-axis coordinate of the coal seam floor min Substitute into the formula: ∑h i =Z max -Z min -M, calculate the total thickness ∑h of the collapsed rock layer. i Let the coefficient of fragmentation be B, and then calculate the accumulation height h at that coordinate. c =B∑h i Collapse height h k = Z max -Z min ; b3. When the roof accumulation height of the goaf is h c Greater than or equal to the collapse height h k When the accumulated gangue reaches a certain height, it indicates that the accumulated gangue has come into contact with the roof and formed an effective support for the roof. In FLAC3D software, the material model along the entire z-direction at this coordinate within the goaf is changed to a Double Yield model, and the basic element is defined as a Fill group. c Less than the collapse height h k At this time, the piled gangue has not yet formed an effective support for the roof and does not provide any force. Therefore, it has no impact on the overall model calculation effect and the goaf is not treated. b4. Repeat steps b1 to b3 above to iterate and judge the coordinates in the entire goaf space, thereby realizing the numerical simulation of the gangue accumulation process during the simulation of mining in the entire geological model.
3. The method for numerical inversion of stress evolution in mining-induced strata and equivalent transformation of stress paths according to claim 1, characterized in that... The steps of the numerical implementation method for the compaction bearing capacity of stockpiled gangue are as follows: c1. Plastic hardening parameters for the Double Yield material model simulating accumulated gangue: coefficient R, cap pressure p c Plastic strain increment ε p Perform calibration; c2. Using the fish language, mechanical parameters are assigned to the basic elements in the Fill group to simulate the compaction bearing mechanical characteristics of gangue in the goaf.
4. The numerical inversion method for stress evolution of mining-induced rock strata and equivalent transformation of stress paths as described in claim 3, characterized in that... Step c1 specifically includes: c11. Sample the crushed gangue piled up in the goaf area in the simulated area, and conduct relevant mechanical property tests on the crushed gangue piled up to obtain the basic mechanical parameters of the crushed gangue piled up: density ρ, cohesion c, internal friction angle φ, tensile strength σt, bulk modulus K, and shear modulus G. c12. Using a uniaxial compression test device for stockpiled crushed gangue, obtain stress and strain data during the compression process of stockpiled crushed gangue; c13. Calculate the coefficient R in the Double Yield material model, which controls the ratio of axial compression to confining pressure; c14. Characterize the compressive bearing characteristics of piled and crushed gangue using the Salamon constitutive relation, and substitute the stress and strain data obtained from the experiment into the cap pressure p. c and plastic strain increment ε p The calculation formula yields the corresponding values.
5. The numerical inversion method for stress evolution of mining-induced strata and equivalent transformation of stress paths according to claim 4, characterized in that... Calculate the coefficient R in step c13 using the following formula: ； Calculate the cap pressure in step c14 using the following formula: ； The plastic strain increment in step c14 is calculated using the following formula: in Indicates axial pressure, λ represents axial strain; K and G are the bulk modulus and shear modulus, respectively.
6. The method for numerical inversion of stress evolution in mining-induced strata and equivalent transformation of stress paths according to claim 1, characterized in that, The stress path equivalent transformation method refers to using the pressure in three directions of a triaxial testing instrument to equivalently replace the stress obtained from the unit volume monitoring in the numerical model, i.e., σ1=σ z σ²=σ y σ3=σ x The axial compression σ1 is equivalent to the vertical stress σ z The confining pressures σ2 and σ3 are equivalent to the horizontal stress σ y σ x The stress evolution law obtained during coal seam mining is equivalently converted into mining-induced stress path, and the actual mining-induced stress environment is restored by triaxial testing instruments.