A coal seam floor hydraulic fracturing fracture topology quantification and permeability improvement and outburst elimination prediction method

Abstract

The application provides a coal seam floor hydraulic fracturing fracture topology quantification and outburst elimination prediction method, and belongs to the technical field of coal mine gas control and coal seam gas development. The application aims at the problems of single description index of existing floor fracturing fracture network, disconnection between the permeability model and the real fracture network, and inability to dynamically predict the outburst elimination effect in advance, reconstructs a three-dimensional model of the coal seam fracture network through true triaxial hydraulic fracturing physical simulation and industrial CT scanning, introduces the complex network theory to calculate the fractal dimension, connectivity, clustering coefficient and small world coefficient, establishes a comprehensive evaluation system for the quality of the fracture network, constructs a permeability dynamic evolution model embedded with the topological characteristics, realizes the space-time prediction of the outburst elimination effect by combining with the multi-field coupling numerical simulation, and establishes a fracturing parameter-fracture topology-outburst elimination standard period quantitative atlas. The application can quantitatively represent the connectivity characteristics of the fracture network, realizes the advanced prediction and parameter optimization of the permeability improvement and outburst elimination effect, and is suitable for the fracturing engineering design of the deep soft and outburst coal seam floor.

Description

Technical Field

[0001] This invention relates to the field of coal mine gas control and coalbed methane development technology, specifically a method for predicting the permeability enhancement and gas outburst reduction effect of hydraulic fracturing on the coal seam floor by combining complex network theory and fracture network topology quantification. Background Technology

[0002] With the gradual depletion of shallow coal resources, the development of deep coal seams has become a strategic direction for ensuring my country's energy security. Deep coal seams are generally located in complex geological environments characterized by high ground stress, high gas pressure, and low permeability, leading to a sharp increase in the risk of dynamic disasters such as coal and gas outbursts. Cross-layer fracturing technology in the coal seam floor strata involves drilling fracturing boreholes in the coal seam floor strata and performing directional or conventional hydraulic fracturing from bottom to top. This drives fractures upwards through the coal-rock interface and into the coal seam, forming a complex network of gas-enhancing fractures within the target coal seam. This technology effectively avoids the difficulties of borehole formation and the risk of borehole collapse in soft coal seams, and has become an important means of overcoming the bottleneck of permeability enhancement in deep, soft coal seams.

[0003] However, the fracture network formed by cross-layer fracturing in the base plate has multi-scale and strong heterogeneity characteristics, and its connectivity efficiency and spatial configuration directly determine the permeability enhancement and outburst reduction effect. The existing technology has the following shortcomings: (1) The existing fracture network description index is singular, relying only on simple geometric parameters such as fracture length, width, and density, and cannot extract topological features such as fractal dimension, connectivity, and clustering coefficient, making it difficult to quantitatively characterize the connectivity efficiency and spatial filling capacity of the permeability enhancement fracture network; (2) The evaluation of fracturing effect lacks a multi-dimensional quantitative index system that integrates fracture scale, connectivity, and distribution density, and the mapping relationship between the base plate fracturing construction parameters and the quality of the fracture network formed in the coal seam is unclear, resulting in a lack of reliable basis for fracturing design optimization; (3) The existing permeability The permeability evolution model is mostly based on the idealized parallel plate fracture assumption, which is seriously disconnected from the real complex fracture network structure formed by bottom plate fracturing. It fails to introduce the topological characteristics of the fracture network as the core control variable and cannot reflect the lasting control effect of the fracture network topology on the flow capacity; (4) The prediction of the gas outburst reduction effect depends on the field measured indicators after the extraction reaches the standard for post-event verification. It is difficult to dynamically predict the spatiotemporal evolution of stress field, gas pressure field, permeability field and effective gas outburst reduction range before or during the extraction. It lacks the ability to conduct forward-looking analysis on the dynamic mechanism of bottom plate fracturing permeability enhancement to promote gas outburst reduction. Therefore, a new method for predicting the gas outburst reduction effect of hydraulic fracturing permeability enhancement in coal seam bottom plate is needed. Summary of the Invention

[0004] To address the aforementioned technical shortcomings, the purpose of this invention is to provide a method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in coal seam floor, which quantitatively characterizes the connectivity of the fracture network and enables advanced prediction and parameter optimization of permeability enhancement and outburst reduction effects.

[0005] To achieve the above objectives, the present invention adopts the following technical solution:

[0006] A method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in coal seam floor includes the following steps:

[0007] Step 1: Prepare coal-rock composite samples, build a true triaxial hydraulic fracturing experimental system integrating loading, pumping and monitoring functions, carry out physical simulation experiments of hydraulic fracturing of coal seam floor, simultaneously collect pumping pressure and acoustic emission signals and track the spatiotemporal evolution trajectory of fracture propagation, perform industrial CT scans on samples before and after fracturing, and reconstruct a three-dimensional digital model of the internal fracture network of the coal seam.

[0008] Step 2: Abstract the three-dimensional fracture network into a complex network graph model, calculate multi-dimensional topology parameters, use the TOPSIS method to determine the weights of each topology parameter, and establish a weighted comprehensive scoring model for fracture network quality.

[0009] Step 3: Set up the fracturing parameter set, design a multi-factor orthogonal experimental scheme, obtain the fracture network topology parameters under different combinations of fracturing parameters, determine the main controlling factors through range analysis and variance analysis, establish a quantitative mapping model between fracturing parameters and fracture network comprehensive score, and complete parameter optimization.

[0010] Step 4: Abstract the fracturing coal seam into a dual-pore medium, construct an initial permeability mapping function with topological parameters as the core control variables, couple the effective stress and deabsorption shrinkage effect to establish a dynamic evolution equation for permeability and complete parameter calibration;

[0011] Step 5: Construct a stress-seepage multi-field coupled mathematical model, embed the fracture network of the optimized fracturing scheme, perform finite element numerical solution, dynamically identify and divide the anti-outburst region based on the simulation results, and output spatiotemporal evolution data of the anti-outburst effect.

[0012] Step six: Conduct parameter sensitivity analysis to identify key control parameters, construct a three-dimensional quantitative relationship map of fracturing parameters, fracture topology, and anti-outburst cycle, and select the optimal engineering application scheme through reverse query.

[0013] Preferably, step one specifically includes:

[0014] Sandstone samples from the bottom of the target coal seam and coal samples from the coal seam were collected, processed into cubic specimens, and cemented into a coal-rock composite. Fracturing boreholes were drilled in the sandstone layer and a simulated wellbore was inserted. Triaxial unequal principal stress was applied to simulate the in-situ in-situ stress state. Fracturing fluid containing fluorescent tracers was injected into the simulated wellbore, and the pumping pressure curve and acoustic emission signal were recorded simultaneously. A time-difference positioning algorithm was used to perform three-dimensional dynamic positioning of micro-fracture events. High-precision industrial CT scans were performed on the specimens before and after fracturing. The spatial coordinates, geometric dimensions, aperture, and distribution orientation of newly formed fractures were extracted through threshold segmentation and three-dimensional reconstruction to establish a three-dimensional digital model of the fracture network.

[0015] Preferably, the multi-dimensional topological parameters in step two include fractal dimension, connectivity, clustering coefficient, and small-world coefficient;

[0016] The fractal dimension is calculated using box counting to quantify the complexity of spatial filling and multi-scale self-similarity in the fracture network; the connectivity rate is the ratio of the number of nodes in the largest interconnected component to the total number of nodes in the fracture network, reflecting the overall connectivity efficiency of the fracture network; the clustering coefficient is the average probability that the neighbors of nodes in the network are neighbors of each other, characterizing the degree of development of closed-loop structures in the fracture network; the small-world coefficient is calculated by comparing the characteristic path length and clustering coefficient of the actual network with those of a random network of the same size, judging the transmission efficiency and local connectivity characteristics of the fracture network.

[0017] Preferably, step three specifically includes:

[0018] The injection flow rate, fracturing fluid viscosity, angle between the branch hole and the coal-rock interface, and branch hole spacing were selected as the fracturing parameter set. A multi-level parameter value was set to design a multi-factor orthogonal experimental scheme. Physical simulation experiments were carried out in sequence according to the orthogonal experimental scheme, and the fracture network topology parameters and comprehensive scores corresponding to each experimental group were calculated. Range analysis and variance analysis were performed on the experimental results to determine the significance ranking and main controlling factors of each fracturing parameter on the fracture network quality. A quantitative mapping model between fracturing parameters and fracture network comprehensive score was established using multivariate nonlinear regression. Parameter optimization was carried out with the goal of maximizing the fracture network comprehensive score, and the optimal fracturing scheme was output.

[0019] Preferably, step four specifically includes:

[0020] The coal seam after floor fracturing is abstracted as a dual-porosity medium composed of a fracture system and a coal matrix block system; the fractal dimension and connectivity are used as core topological control variables, and an initial permeability mapping function for the fracture network is constructed by combining fracture porosity and maximum fracture length:

[0021] In the formula, The fissure porosity is obtained through volumetric statistics of CT images in step one. To determine the maximum crack length, extract the lengths of all crack segments from step one and take the maximum value. The fitting coefficient is set to 0.2.

[0022] By introducing the dynamic disturbance of fracture aperture caused by effective stress change and gas desorption shrinkage, a dynamic evolution equation for permeability is established:

[0023] In the formula, The fracture compressibility coefficient; Effective stress; The initial effective stress before the start of extraction; The matrix shrinkage coefficient; To solve the absorption of contraction strain; Described by Langmuir-type equations;

[0024] Using measured permeability data under different fracturing conditions, the model parameters were inverted and calibrated, and the accuracy was verified.

[0025] Preferably, step five specifically includes:

[0026] By combining the control equations of coal and rock deformation field, gas seepage field, and gas desorption and diffusion equation, a multi-field coupled mathematical model of stress and seepage is constructed for the entire process of coal seam floor fracturing and extraction. Based on the actual geological conditions of the mine, an engineering-scale geometric model of coal seam-floor strata is established, and the fracture network of the optimized fracturing scheme is embedded into the geometric model. The ground stress boundary, gas pressure boundary, and extraction borehole negative pressure boundary are set. The multi-field coupled model is solved using the finite element method, and the spatiotemporal distribution data of coal seam stress field, gas pressure field, and permeability field during the extraction process are output. Based on the critical value of outburst risk, the outburst reduction zone is dynamically identified and divided, and the outburst reduction range and the time to achieve the outburst reduction target are output for different extraction times.

[0027] Preferably, step six specifically includes:

[0028] By changing individual fracturing and extraction parameters while keeping other parameters constant, the sensitivity of each parameter to the impact on the peak burst reduction period and range is analyzed by comparing simulation results, and key control parameters are identified. Parameter scanning simulations are then performed on the key control parameters within their value ranges, and the simulation results are correlated with the input parameters to create a three-dimensional quantitative relationship map of fracturing parameters, fracture topology, and peak burst reduction period. Based on the engineering requirements for the peak burst reduction range and period, the optimal fracturing process parameters and extraction system layout are selected through reverse lookup using the quantitative map, forming an application scheme for bottom plate fracturing permeability enhancement and peak burst reduction engineering technology.

[0029] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0030] 1. This invention proposes a multi-index quantitative characterization system for the topology of fracture networks by introducing complex network theory, with fractal dimension, connectivity, clustering coefficient and small-world effect coefficient as the core.

[0031] 2. This invention establishes a multi-index weighted comprehensive scoring model by integrating fracture scale characteristics, connectivity performance, and transmission efficiency.

[0032] 3. This invention establishes a quantitative mapping model between fracturing construction parameters and fracture network topology parameters through systematic physical simulation and regression analysis. The optimal construction parameters can be deduced based on the target fracture network quality, significantly shortening the parameter optimization cycle.

[0033] 4. This invention constructs a closed-loop coupled model of fracture topology-permeability evolution by embedding fractal dimension and connectivity as core independent variables into the initial permeability mapping function and coupling the dynamic influence of effective stress change and matrix deabsorption shrinkage.

[0034] 5. This invention establishes a positive prediction method for the spatiotemporal evolution of stress field, gas pressure field, permeability field, and gas outburst elimination zone throughout the entire extraction cycle by deeply integrating fracture network topology characterization, permeability dynamic evolution model, and multi-field coupled numerical simulation.

[0035] 6. Through extensive parameter sensitivity analysis, this invention constructs a quantitative relationship map of fracturing parameters, fracture topology, and burst elimination cycle, enabling reverse querying of engineering objectives and coordinated optimization of fracturing and extraction. Attached Figure Description

[0036] Figure 1 This is a flowchart illustrating the overall technical process of the method of the present invention.

[0037] Figure 2 This is a schematic diagram showing the distribution of hydraulic fracturing fractures in the coal seam floor.

[0038] Figure 3 A schematic diagram illustrating the abstraction of a fracture network into a complex network graph model.

[0039] in:

[0040] 1. Coal seam; 2. Coal seam fracture; 3. Coal-rock interface; 4. Rock stratum fracture; 5. Rock stratum; 6. Fracturing hole; 7. Actual fracture; 8. Actual fracture intersection; 9. Edge of complex network model; 10. Node of complex network model. Detailed Implementation

[0041] The invention will now be further described with reference to the accompanying drawings.

[0042] like Figure 1 As shown, a method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor includes the following steps:

[0043] Step 1: Physical Simulation of Coal Seam Floor Fracturing and 3D Reconstruction of Fracturation Network

[0044] like Figure 2 As shown, 1 represents the coal seam; 2 represents the coal seam fracture; 3 represents the coal-rock interface; 4 represents the rock strata fracture; 5 represents the rock strata; and 6 represents the hydraulic fracturing hole.

[0045] 1. Preparation of coal-rock composite samples

[0046] Sandstone samples from the bottom of the target coal seam and coal samples from the soft, outburst-prone coal seam were collected from the site. Both materials were processed into 150mm × 150mm × 150mm cubes. High-strength epoxy resin was used to bond the sandstone blocks and coal blocks to form a coal-rock composite sample, ensuring a smooth coal-rock interface and controllable bonding strength. A borehole (10mm diameter, 10-15mm depth from the interface) was drilled in the center of the sandstone layer along a direction perpendicular to the interface, simulating a floor fracturing borehole. A simulated wellbore with an outer diameter of 10mm was inserted, and epoxy resin was used to seal the gap between the wellbore and the borehole wall.

[0047] 2. True Triaxial Hydraulic Fracturing Physical Simulation Experiment

[0048] A small-scale hydraulic fracturing experimental system for coal-rock assemblies was constructed, integrating a true triaxial loading and servo control subsystem, a hydraulic fracturing pump injection subsystem, and a real-time acoustic emission monitoring and positioning subsystem. The sample was placed in a true triaxial loading chamber, and triaxial unequal principal stresses were applied to simulate in-situ geostress. Fracturing fluid containing a fluorescent tracer was injected into the simulated wellbore via a high-pressure, constant-speed, constant-pressure pump, while the pump injection pressure curve and acoustic emission signal were simultaneously monitored and recorded in real time.

[0049] This invention employs the following steps to refine the original AE signal in order to obtain the spatiotemporal evolution law of crack propagation:

[0050] (1) Signal preprocessing and filtering

[0051] The raw AE signal contains environmental noise. First, a bandpass filter is applied to remove low-frequency mechanical vibrations and high-frequency electrical noise. Then, short-time averaging and long-time averaging algorithms are used to identify valid AE events. The STA window length is set to 0.1ms (100 sampling points), and the LTA window length is set to 1ms (1000 sampling points). When the STA / LTA ratio exceeds the set trigger threshold, it is determined as the start of a valid micro-fracture event.

[0052] (2) Accurate pickup of P wave arrival time

[0053] For each triggered AE event, the arrival time of the P wave in each channel is automatically picked up using the Akaike Information Criterion (AIC).

[0054] (3) Localization of three-dimensional micro-fracture sources

[0055] Based on the known three-dimensional dimensions of the sample and the spatial coordinates of the sensor and the arrival time of the first P wave. A time-difference positioning equation set was constructed. Assuming the coal-rock assemblage is a homogeneous equivalent body with indistinct anisotropy, the longitudinal wave velocity was taken as... The Geiger iterative method is used for nonlinear least squares solution, with the objective function being to minimize the sum of squared location residuals.

[0056] In the formula, The coordinates of the micro-fracture source to be located are shown. The moment of tremor is defined. The iteration termination condition is set to a position update amount less than 0.1 mm or a residual change rate less than 1e-6. The final output is the three-dimensional spatial coordinates of each AE event. And when it occurs.

[0057] (4) Construction of crack evolution trajectory and spatiotemporal evolution diagram

[0058] Arrange all the located event points in chronological order to form a four-dimensional point cloud. The DBSCAN clustering algorithm was used to remove discrete noise points. For each cluster, the convex hull algorithm was used to extract the contour of the fracture propagation front according to the time series, and a fracture spatiotemporal evolution diagram was drawn to distinguish the three key stages of floor initiation, interface crossing, and coal seam network formation.

[0059] 3. CT scans and three-dimensional reconstruction of fracture networks before and after fracturing

[0060] Before and after the fracturing experiment, the samples were subjected to high-precision industrial CT scans. Based on the CT images before and after fracturing, the spatial coordinates, geometric dimensions, aperture, and spatial distribution of newly formed fractures were extracted using threshold segmentation and 3D reconstruction techniques. This generated a 3D digital model of the fracture network within the coal seam, providing a data foundation for subsequent topology quantization.

[0061] (1) Image preprocessing: The original CT slices were denoised using a nonlocal mean filtering algorithm to remove noise while preserving the edge features of the crack.

[0062] (2) The grayscale image is divided into four categories using the multi-threshold OTSU algorithm: air / fracture (lowest grayscale), coal matrix (lower grayscale), sandstone matrix (medium grayscale), and high-density minerals (highest grayscale). The lowest grayscale class is marked as "fracture phase" to generate an initial binary image.

[0063] (3) Morphological opening operation is used to remove isolated noise points (structural element: 3×3 pixel disk), morphological closing operation is used to repair the interruption of the fracture connectivity, and isolated regions with a volume of less than 50 voxels are deleted by the connected component labeling algorithm to obtain a clean fracture network binary image.

[0064] (4) Stack all processed binary slices in spatial order and use the moving cube algorithm to generate a three-dimensional surface model of the triangular facets of the fracture network. On this basis, use the three-dimensional topology refinement algorithm to skeletonize the fracture network volume data and output a three-dimensional skeleton line model containing the center lines of all fracture segments for subsequent topology quantization data input.

[0065] Step 2: Quantification of Slit Network Topology Features and Construction of a Multi-Indicator Evaluation System

[0066] 1. The fracture network is abstracted into a complex network graph model.

[0067] like Figure 3 As shown, 7 represents the actual crack; 8 represents the intersection of the actual crack; 9 represents the edge of the complex network model; and 10 represents the node of the complex network model.

[0068] The 3D fracture network skeleton line data output from step one is abstracted into an undirected weighted graph:

[0069] In the formula, A set of nodes; Let it be the set of edges; Let be the set of edge weights. In this embodiment, the length is taken as the weight.

[0070] The nodes and edges are defined as follows:

[0071] Connecting Node: In the skeleton line model, if three or more crack segments intersect at a point, that point is defined as a connecting node.

[0072] Endpoint: An endpoint is defined as one end of a crack segment that is connected to only one edge.

[0073] Isolated nodes: If a voxel belongs to a fissure but is not connected to any other fissure voxels, it will be ignored or marked separately in subsequent connectivity calculations.

[0074] A continuous crack segment between any two adjacent nodes is defined as an edge. The properties of each edge include:

[0075] length Euclidean distance (mm) between the two endpoints.

[0076] Average opening : The average local opening (mm) of all skeleton points along this edge.

[0077] Azimuth : The projection direction angle of the edge onto the horizontal plane (0~360°).

[0078] inclination : The angle between the side and the horizontal plane (-90°~90°).

[0079] This embodiment uses a complex network analysis tool to construct a graph model and store the attribute data of nodes and edges.

[0080] 2. Calculation of topology parameters

[0081] Based on the constructed graph model, multi-dimensional topological parameters are calculated to characterize the connectivity efficiency and space-filling capacity of the gap network, specifically including:

[0082] (1) Fractal dimension

[0083] The fractal dimension is calculated using box counting. By covering the fracture network with cubic boxes of different sizes, the number of non-empty boxes required to cover the skeleton is counted, and the fractal dimension is obtained by log-linear regression slope. This method is used to quantify the complexity of the fracture network's space-filling and its multi-scale self-similarity.

[0084] Set the box side length The value range is: 2, 4, 8, 16, 32, 64, 128 pixels. The entire sample space is divided into sections with side lengths of... The cubic mesh is used. Each mesh is traversed, and it is determined whether it contains at least one fracture skeleton voxel (read from the binarized fracture volume data in step one). The number of non-empty meshes is counted. .

[0085] Data fitting: calculation , A least-squares linear fit is performed on the data points. The slope is the fractal dimension. Since the crack network in this study is a three-dimensional volume, The value is between 2 and 3. To facilitate the construction of subsequent permeability models, the fractal dimension is defined as the volume-filling dimension, with a value ranging from 2 to 3. The larger the value, the more abundant the small-scale fractures.

[0086] (2) Connectivity

[0087] Connectivity is defined as the ratio of the number of nodes with the largest number of interconnections in a split network to the total number of nodes in the network. It is used to reflect the overall connectivity efficiency of a split network.

[0088] Use complex network analysis tools to obtain a list of nodes for all connected components. Find the largest connected component and denote its number of nodes as . Total number of nodes Includes all nodes. Calculate connectivity: If the fracture network is completely isolated, then If fully connected, then .

[0089] (3) Clustering coefficient

[0090] Clustering coefficients reflect the degree of development of triangular closed-loop structures in a fracture network by calculating the average probability that the neighbors of nodes in the network are also neighbors of each other.

[0091] Local clustering coefficient of a node: for a node (Spend The actual number of edges between its neighboring nodes is .but:

[0092] like ,but .

[0093] Average clustering coefficient of the network:

[0094] in This represents the total number of nodes in the network.

[0095] The average clustering coefficient of the network can be obtained by calculating the clustering coefficient of each node using complex network analysis tools and then averaging the results. .

[0096] (4) Small-world effect coefficient

[0097] Small-world networks exhibit both high clustering coefficients and short characteristic path lengths. By combining characteristic path length and clustering coefficients, a method can be used to determine whether a gap network possesses both efficient transmission and locally dense connectivity. Characteristic path length is defined as the average of the shortest path lengths between all node pairs. For disconnected graphs, typically only node pairs within the most connected component are calculated.

[0098] Calculate the characteristic path length of a real network and clustering coefficient Generate at least 100 nodes of equal size (N nodes, edge connection probability is 1). An Erdős–Rényi random network instance was generated, and the average characteristic path length of the random network was calculated. and average random network clustering coefficient .

[0099] Small World Coefficient :

[0100] like If the value is large enough, the gap network exhibits a small-world effect; a larger value indicates a more pronounced combination of efficient transmission and locally dense connectivity. Generally, For the network of the weak world, For strong small-world networks. In this study The typical value is 1.2~2.5.

[0101] 3. Multi-indicator comprehensive scoring model

[0102] This invention uses the TOPSIS method to determine the weights of fractal dimension, connectivity, clustering coefficient, and small-world effect coefficient, and establishes a weighted comprehensive scoring model to quantitatively and comprehensively evaluate the quality of the fracture network formed by coal seam floor fracturing.

[0103] (1) There is One fracturing test specimen to be evaluated. Several evaluation metrics (fractal dimension, connectivity, clustering coefficient, small-world effect coefficient). The original data matrix is ​​constructed as follows: ,in Indicates the first The first scheme Individual indicator values.

[0104] (2) Since the four topological parameters have different dimensions and numerical ranges, they need to be normalized to the [0,1] interval. This invention uses vector normalization:

[0105] In the formula, ; Obtain the normalized matrix .

[0106] (3) Based on the nine sets of data from the orthogonal experiment in step three, the entropy weight method was used to calculate the weights of each parameter, and the weight vector was obtained. Calculate the weighted normalized value. The weighted matrix is ​​obtained. .

[0107] (4) Determine the ideal solution and negative ideal solution Since all indicators are benefit-oriented, the calculations are as follows: , ,

[0108] (5) Calculate the Euclidean distance from each scheme to the positive and negative ideal solutions.

[0109] Distance to the ideal solution:

[0110] Distance to the negative ideal solution:

[0111] (6) Calculate the relative closeness

[0112] The larger the value, the higher the value. The better the solution. Sort all fracturing schemes by value from largest to smallest, and select... The largest possible solution is considered the optimal fracturing solution under the current conditions.

[0113] All four indicators mentioned above are benefit-oriented indicators: the higher the fractal dimension, the greater the connectivity, the higher the clustering coefficient, and the more significant the small-world effect, the better the quality of the fracture network.

[0114] Step 3: Analysis of the impact of fracturing parameters on fracture network quality and optimization of the scheme

[0115] 1. Fracturing parameter set design

[0116] Based on engineering practice and preliminary sensitivity analysis of cross-seam fracturing in the coal seam floor, the following four key fracturing parameters were selected as influencing factors: injection flow rate. (mL / min), fracturing fluid viscosity (mPa·s), angle between the branch hole and the coal-rock interface (°), branch hole spacing (mm).

[0117] 2. Multi-factor orthogonal experimental design

[0118] This embodiment employs an orthogonal design of main effects and does not specifically examine interactions. In practical applications, local interaction experiments can be conducted for verification after preliminary optimization. A three-level design is used to examine potential nonlinear effects. The level values ​​are determined based on engineering experience.

[0119] Level 1: Injection flow rate 3 mL / min; fracturing fluid viscosity 1 mPa·s; branch hole angle 90°; branch hole spacing 30 mm. Level 2: Injection flow rate 6 mL / min; fracturing fluid viscosity 10 mPa·s; branch hole angle 60°; branch hole spacing 40 mm. Level 3: Injection flow rate 12 mL / min; fracturing fluid viscosity 100 mPa·s; branch hole angle 30°; branch hole spacing 50 mm.

[0120] Nine orthogonal experiments were conducted, with injection flow rate as the primary variable. True triaxial hydraulic fracturing physical simulation experiments were performed strictly according to the method described in step one, and the comprehensive quality score of the fracture network for each group of experiments was calculated based on the comprehensive scoring method described in step two. ( Simultaneously, the physical simulation experiments described in step one are carried out sequentially according to the plan to obtain the fracture network topology parameters under different combinations of fracturing parameters.

[0121] 3. Regression Analysis and Mapping Model

[0122] Range analysis was performed on the orthogonal experimental results to calculate the range values ​​corresponding to each factor. :

[0123] In the formula, The maximum average of the factor scores. This represents the minimum average score of the factor. A larger range indicates a greater impact of that factor on the overall quality score of the fracture network. The more significant the impact.

[0124] Perform analysis of variance on the results of the orthogonal experiment and calculate the total sum of squares. :

[0125] In the formula, This represents the overall mean of the factor scores.

[0126] Calculate the sum of squared deviations for each factor. :

[0127] In the formula, For horizontal numbers. Sum of squared deviations of errors. :

[0128] In the formula, Injection flow rate Sum of squared deviations fracturing fluid viscosity Sum of squared deviations The angle between the branch hole and the coal-rock interface Sum of squared deviations Branch hole spacing Sum of squared deviations.

[0129] Calculate the F-value:

[0130] The significance of each factor's influence was quantitatively determined using the F-test to identify the controlling factor. Finally, multiple nonlinear regression was employed to establish a quantitative mapping model between fracturing parameters and the comprehensive fracture network score.

[0131] Step 4: Constructing a dynamic evolution model of penetration rate embedded with the topological features of the fracture network

[0132] 1. Dual-porosity medium hypothesis

[0133] The coal seam after the bottom plate is fractured is regarded as a dual porous medium consisting of a fracture system and a coal matrix block system, in which the fracture system is the main seepage channel and the coal matrix block is the main gas storage space.

[0134] The geometric model makes the following assumptions: the fracture system can be equivalent to a set of mutually orthogonal parallel plate fracture networks; the coal matrix block has a side length of... The cube, the gap spacing ; fracture aperture Much smaller than the fracture spacing. The irregularity, connectivity efficiency, and space-filling capacity of the actual fracture network are determined by the fractal dimension calculated in step two. and connectivity Make corrections.

[0135] 2. Mapping model between initial permeability and topological parameters

[0136] This model is modified based on the permeability theory of fractal porous media and the volume-filling characteristics of a three-dimensional fracture network. The total permeability of the fracture network can be expressed as the sum of contributions from fractures at various scales. The fractal dimension calculated in step two... With connectivity As a core control variable, the initial penetration rate of the fracture network is constructed. Mapping function:

[0137] In the formula, The fissure porosity is obtained through volumetric statistics of CT images in step one. To determine the maximum crack length, extract the lengths of all crack segments from the skeleton data in step one and take the maximum value. The fitting coefficient is set to 0.2. This function reflects the persistent control effect of the mesh topology on the initial flow conduction capability.

[0138] The coefficient 3 in the index comes from the relationship between permeability and volumetric strain:

[0139] In the formula, The value represents volumetric strain. Both effective stress variation and deabsorption shrinkage affect permeability through volumetric strain; therefore, a coefficient of 3 is uniformly adopted. This form is consistent with the Palmer & Mansoori model.

[0140] 3. Dynamic evolution equation of permeability during extraction process

[0141] During gas extraction, the dynamic change in permeability is controlled by the effective stress effect and the gas desorption shrinkage effect: the increase of effective stress leads to the closure of fractures. As the gas pressure decreases, the effective stress acting on the coal and rock skeleton increases, the fractures are compressed, the opening decreases, and the permeability decreases. Gas desorption causes matrix shrinkage. After gas desorption, the coal matrix shrinks, the fracture opening increases, and the permeability increases.

[0142] (1) Effective stress effect

[0143] According to Terzaghi's effective stress principle:

[0144] In the formula: The total stress is the sum of the pressure of the overlying strata and the tectonic stress; This is the Biot coefficient, set to 0.8; This refers to the pore gas pressure.

[0145] The response of effective stress to crack aperture is usually assumed to be exponential:

[0146] In the formula: This represents the initial crack aperture; The fracture compressibility factor is taken as 0.3 MPa⁻¹; This represents the initial effective stress before extraction begins.

[0147] Since permeability is proportional to the square of aperture, we have:

[0148] (2) Gas desorption and shrinkage effect

[0149] Coal matrix expands after adsorbing methane gas and contracts upon desorption. Matrix shrinkage strain. It can be described by Langmuir-type equations:

[0150] In the formula: The maximum shrinkage strain of Langmuir is taken as 0.01; is the Langmuir pressure constant, which corresponds to the gas pressure when the adsorption capacity reaches 0.5 times the maximum adsorption capacity, and is taken as 1 MPa.

[0151] An increase in shrinkage strain leads to an increase in fracture aperture. Assuming a linear relationship between fracture aperture and matrix shrinkage strain, with a coefficient of... :

[0152] In the formula, The relative change in crack aperture caused by unit shrinkage strain is taken as 6.

[0153] From the relationship between permeability and the square of aperture, we get:

[0154] (3) Coupled dynamic evolution equations

[0155] Considering both effective stress and deabsorption shrinkage, the change in total permeability is the product of the two:

[0156] To simplify the coefficients in the exponent, the commonly used form is:

[0157] The coefficient 3 is used here to align with the volumetric strain-permeability relationship (an exponential relationship between permeability and volumetric strain) found in common literature. Differences between coefficients 2 and 3 can be absorbed during the calibration process. This invention uses coefficient 3 to maintain compatibility with most models.

[0158] 4. Experimental calibration of model parameters

[0159] Using the measured permeability data under different fracturing conditions in step one, the unknown parameters in the above model are inverted and calibrated to ensure the prediction accuracy of the model.

[0160] (1) Initial permeability mapping model parameters calibration

[0161] Using at least nine sets of samples under different fracturing parameters obtained from the orthogonal experiments in step three, the fracture porosity of each sample was obtained through step one. fractal dimension Connectivity and maximum crack length The initial permeability of each sample was obtained by steady-state permeability testing. .

[0162] Construct the dataset and fit the formula:

[0163] Using the least squares method, let:

[0164] get :

[0165] (2) Fracture compressibility coefficient calibration

[0166] Using permeability data with varying effective stress (fixed gas pressure, changed confining pressure), an exponential relationship was fitted:

[0167] Take the logarithm of both sides:

[0168] The slope of the linear regression is The result can be obtained. .

[0169] For the desorption-permeability effect, a desorption-permeability experiment was conducted separately: While maintaining a constant confining pressure, the gas pressure was gradually decreased, and the steady-state permeability at each pressure was recorded. The relationship between permeability and gas pressure was then fitted using a model.

[0170] in The only unknown parameter was obtained from an independent isothermal adsorption-expansion experiment. Calibration can be achieved through nonlinear regression:

[0171] This can be obtained through nonlinear fitting. The typical value range is 2 to 10.

[0172] (3) Langmuir shrinkage strain parameters and Measurement

[0173] Standard samples were processed from unfractured raw coal samples and placed in an isothermal adsorption-expansion synchronous measurement device. Methane gas at different pressures was introduced, and after equilibrium was reached, the adsorption volume and volumetric strain were recorded simultaneously.

[0174] The adsorption data were fitted using the Langmuir equation:

[0175] The strain data were fitted using Langmuir-type strain equations:

[0176] Determined step by step through nonlinear regression (Levenberg-Marquardt algorithm) , and .

[0177] (4) Model validation

[0178] Leave-one-out cross-validation was employed: each time, one set of orthogonal experimental data was retained as the validation set, and the model parameters were recalibrated using the remaining data to predict the penetration rate of that set, and the relative error was calculated. If the average relative error of the nine validations was less than 20%, the model calibration accuracy was considered to meet engineering requirements. In addition, an extra set of orthogonal experimental parameter combinations not involved in the regression was performed, and the model predictions were compared with the measured values ​​to further verify the model's generalization ability.

[0179] Step 5: Multi-field coupled numerical simulation and spatiotemporal prediction of conflict mitigation effects throughout the extraction process

[0180] 1. Multi-field coupled control equations

[0181] Based on the permeability dynamic evolution model established in step four, and combined with the coal and rock deformation field control equation, the gas seepage field control equation, and the gas desorption and diffusion equation, a "stress-seepage" multi-field coupled mathematical model for the entire process of coal seam floor fracturing and extraction is constructed.

[0182] This model comprises three physical fields: the coal and rock deformation field, the gas seepage field, and the gas diffusion field. These three fields are coupled together through effective stress, adsorption expansion strain, and crossflow source terms.

[0183] (1) Control equations of coal and rock deformation field

[0184] Based on Biot's theory of porosity elasticity, and considering the expansion strain caused by gas adsorption, the equilibrium equation for the coal-rock mass is:

[0185] The constitutive equation of the coal-rock mass is:

[0186]

[0187]

[0188]

[0189] In the formula: This is the total stress tensor; For strain tensor; For displacement components; Shear modulus; Let Lamé's constant be denoted by . It is the elastic modulus; Poisson's ratio; This is the Biot coefficient, ranging from 0.6 to 0.9; The gas pressure in the fissure; Bulk modulus; For adsorption expansion strain, where For the maximum expansion strain, is the Langmuir pressure constant.

[0190] (2) Control equations of gas seepage field

[0191] The flow of gas in a fractured system follows Darcy's law, and the mass conservation equation is:

[0192] In the formula: The permeability is the fracture permeability. The dynamic viscosity of gas; Gas density; This refers to the porosity of the fractured area. This represents the gas desorption-diffusion source term from the matrix to the fractures in a unit volume of coal and rock.

[0193] (3) Gas desorption-diffusion equation

[0194] Gas desorbed from coalbed methane blocks diffuses into fractures. A quasi-steady-state diffusion model is adopted, with the source term being:

[0195] In the formula: The adsorption time constant; This represents the current average gas content in the matrix; For the volume of Langmuir; This is a factor that converts volume content into mass density.

[0196] In this patent, for ease of numerical implementation, the equivalent cross-flow coefficient method is recommended, which equates matrix desorption and diffusion to a source term proportional to the pressure difference:

[0197] In the formula The desorption time constant is This represents the average gas content in the matrix.

[0198] For simplicity, desorption can be considered as an instantaneous equilibrium, that is, assuming that the gas content in the matrix is ​​always in equilibrium with the fracture pressure. It can be expressed as a function of the rate of change of pressure. However, in order to improve the prediction accuracy, this invention adopts a dual porosity model to solve for the fracture pressure and matrix content separately.

[0199] (4) Summary of Coupling Relationships

[0200] The deformation field affects the fracture aperture through effective stress and adsorption strain, thereby altering permeability. ;

[0201] Changes in fissure pressure cause changes in effective stress, which in turn affect the deformation field;

[0202] Channeling between matrix and fracture It depends on the difference between fracture pressure and matrix content;

[0203] Deabsorption and contraction strain Coupled to the deformation field.

[0204] 2. Simulate setting conditions

[0205] Based on the actual geological conditions and mining conditions of the mine, an engineering-scale geometric model of the coal seam-floor strata combination (200m×100m×40m) was established. The fracture network distribution results generated by the optimized fracturing scheme obtained in step three were then incorporated into the model using an equivalent continuous medium or discrete fracture network method.

[0206] Stress boundary: The top surface is subjected to an equivalent pressure of 20 MPa, the bottom surface is subjected to normal constraint, the four sides are subjected to a pressure of 22 MPa in the X direction and a pressure of 14 MPa in the Y direction.

[0207] Seepage boundary: The top, bottom and all sides are set as non-flow boundaries, and the borehole inner wall is set as a constant pressure boundary. The initial original gas pressure is 2.5 MPa.

[0208] The COMSOL Multiphysics software was used to solve the multi-field coupled model using the finite element method. The spatiotemporal distribution data of the coal seam stress field, gas pressure field, and permeability field at different times during the extraction process were simulated and calculated.

[0209] (4) Quantitative identification of dynamic anti-protrusion zone

[0210] Based on the spatiotemporal distribution data of gas pressure and residual gas content output from the simulation, and according to the critical values ​​of outburst risk (such as gas pressure) in the "Detailed Rules for the Prevention and Control of Coal and Gas Outbursts". Residual gas content The simulation area is dynamically identified and finely divided into non-protrusion zones.

[0211] Residual gas content (Unit: m³ / t) Calculated by the following formula:

[0212] In the formula, the first term is free methane, and the second term is adsorbed methane. The apparent density of coal; Gas molar mass; It is the ideal gas constant; This refers to the absolute temperature of the coal seam.

[0213] Output the range of gas depletion, the morphology of the gas depletion zone, and the time required to achieve gas depletion target for different extraction times, and draw a spatiotemporal evolution map of the gas depletion effect.

[0214] Step Six: Establishment of Fracturing-Extraction Co-optimization Strategy and Quantitative Mapping

[0215] 1. Parameter sensitivity analysis

[0216] In the numerical model of step five, individual fracturing parameters (such as injection flow rate and branch hole density) or extraction parameters (such as extraction negative pressure and borehole spacing) are changed respectively, while keeping other parameters unchanged. By comparing the simulation results, the sensitivity of each parameter to the impact of the fracturing target period and the fracturing range is analyzed, and the key control parameters are identified.

[0217] Define sensitivity index The ratio of the relative rate of change of time to achieve target resolution to the relative rate of change of parameters:

[0218] In the formula, As the baseline parameter value, This is the changed parameter value. If Increasing the parameters leads to a longer de-suppression time; Increasing parameters shortens the muzzle clearing time. Larger parameters have a greater impact on the effectiveness of surge suppression. Test results show that the branch hole spacing and injection flow rate are the most sensitive parameters affecting the time to achieve surge suppression targets, followed by the borehole spacing. Therefore, subsequent quantitative mapping should focus on constructing these two parameters.

[0219] 2. Construction of Quantitative Relationship Map

[0220] For the key parameters identified in the sensitivity analysis, injection flow rate and branch hole spacing were selected as variables on the horizontal and vertical axes, respectively. The output variables were: fracture network comprehensive score and time to achieve fracturing target. Within feasible limits, two key parameters were selected for parameter scanning simulation. The simulation results were correlated with the input parameters to create a three-dimensional quantitative relationship map of "fracturing parameters --- fracture topology --- fracturing target achievement period".

[0221] (1) Fracturing parameters - fracture topology contour map

[0222] With the injection flow rate as the X-axis and the branch hole spacing as the Y-axis, the fracture network comprehensive score is calculated. Using the Z value, a filled contour map is plotted, and the contour values ​​are labeled. This map visually demonstrates the fracture network quality achievable by different fracturing schemes.

[0223] (2) Response surface diagram of fracturing parameters - anti-outburst target period

[0224] Similarly, with the injection flow rate as the X-axis and the branch hole spacing as the Y-axis, the time to achieve the target for bubble elimination is... Draw a contour plot for the Z value. Mark it with a solid red line. The contour lines.

[0225] (3) Combination diagram

[0226] Overlay the two images above: the background is... Colored contour lines, superimposed The contour lines (black dashed lines) allow for simultaneous reading on a single graph: for a given injection flow rate and branch hole spacing, both the fracture quality and the time required for fracture suppression can be determined.

[0227] 3. Engineering Scheme Reverse Selection and Optimization

[0228] Based on the actual engineering requirements for the scope and cycle of outburst suppression, a reverse lookup was performed using an established quantitative relationship graph to select the optimal fracturing process parameters and extraction system layout scheme that meet the engineering objectives. Combining the constraints of on-site geological structure and mining technology, a set of engineering technology application schemes for floor fracturing, permeability enhancement, and outburst suppression in soft, outburst-prone coal seams was developed.

Claims

1. A method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor, characterized in that, Includes the following steps: Step 1: Prepare coal-rock composite samples, build a true triaxial hydraulic fracturing experimental system integrating loading, pumping and monitoring functions, carry out physical simulation experiments of hydraulic fracturing of coal seam floor, simultaneously collect pumping pressure and acoustic emission signals and track the spatiotemporal evolution trajectory of fracture propagation, perform industrial CT scans on samples before and after fracturing, and reconstruct a three-dimensional digital model of the internal fracture network of the coal seam. Step 2: Abstract the three-dimensional fracture network into a complex network graph model, calculate multi-dimensional topology parameters, use the TOPSIS method to determine the weights of each topology parameter, and establish a weighted comprehensive scoring model for fracture network quality. Step 3: Set up the fracturing parameter set, design a multi-factor orthogonal experimental scheme, obtain the fracture network topology parameters under different combinations of fracturing parameters, determine the main controlling factors through range analysis and variance analysis, establish a quantitative mapping model between fracturing parameters and fracture network comprehensive score, and complete parameter optimization. Step 4: Abstract the fracturing coal seam into a dual-pore medium, construct an initial permeability mapping function with topological parameters as the core control variables, couple the effective stress and deabsorption shrinkage effect to establish a dynamic evolution equation for permeability and complete parameter calibration; Step 5: Construct a stress-seepage multi-field coupled mathematical model, embed the fracture network of the optimized fracturing scheme, perform finite element numerical solution, dynamically identify and divide the anti-outburst region based on the simulation results, and output spatiotemporal evolution data of the anti-outburst effect. Step six: Conduct parameter sensitivity analysis to identify key control parameters, construct a three-dimensional quantitative relationship map of fracturing parameters, fracture topology, and anti-outburst cycle, and select the optimal engineering application scheme through reverse query.

2. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, Step one specifically includes: Sandstone samples from the bottom of the target coal seam and coal samples from the coal seam were collected, processed into cubic specimens, and cemented into a coal-rock composite. Fracturing boreholes were drilled in the sandstone layer and a simulated wellbore was inserted. Triaxial unequal principal stress was applied to simulate the in-situ in-situ stress state. Fracturing fluid containing fluorescent tracers was injected into the simulated wellbore, and the pumping pressure curve and acoustic emission signal were recorded simultaneously. A time-difference positioning algorithm was used to perform three-dimensional dynamic positioning of micro-fracture events. High-precision industrial CT scans were performed on the specimens before and after fracturing. The spatial coordinates, geometric dimensions, aperture, and distribution orientation of newly formed fractures were extracted through threshold segmentation and three-dimensional reconstruction to establish a three-dimensional digital model of the fracture network.

3. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, The multidimensional topological parameters in step two include fractal dimension, connectivity, clustering coefficient, and small-world coefficient. The fractal dimension is calculated using box counting to quantify the complexity of spatial filling and multi-scale self-similarity in the fracture network; the connectivity rate is the ratio of the number of nodes in the largest interconnected component to the total number of nodes in the fracture network, reflecting the overall connectivity efficiency of the fracture network; the clustering coefficient is the average probability that the neighbors of nodes in the network are neighbors of each other, characterizing the degree of development of closed-loop structures in the fracture network; the small-world coefficient is calculated by comparing the characteristic path length and clustering coefficient of the actual network with those of a random network of the same size, judging the transmission efficiency and local connectivity characteristics of the fracture network.

4. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, Step three specifically includes: The injection flow rate, fracturing fluid viscosity, angle between the branch hole and the coal-rock interface, and branch hole spacing were selected as the fracturing parameter set. A multi-level parameter value was set to design a multi-factor orthogonal experimental scheme. Physical simulation experiments were carried out in sequence according to the orthogonal experimental scheme, and the fracture network topology parameters and comprehensive scores corresponding to each experimental group were calculated. Range analysis and variance analysis were performed on the experimental results to determine the significance ranking and main controlling factors of each fracturing parameter on the fracture network quality. A quantitative mapping model between fracturing parameters and fracture network comprehensive score was established using multivariate nonlinear regression. Parameter optimization was carried out with the goal of maximizing the fracture network comprehensive score, and the optimal fracturing scheme was output.

5. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, Step four specifically includes: The coal seam after floor fracturing is abstracted as a dual-porosity medium composed of a fracture system and a coal matrix block system; the fractal dimension and connectivity are used as core topological control variables, and an initial permeability mapping function for the fracture network is constructed by combining fracture porosity and maximum fracture length: In the formula, The fissure porosity is obtained through volumetric statistics of CT images in step one. To determine the maximum crack length, extract the lengths of all crack segments from step one and take the maximum value. The fitting coefficient is set to 0.

2. By introducing the dynamic disturbance of fracture aperture caused by effective stress change and gas desorption shrinkage, a dynamic evolution equation for permeability is established: In the formula, The fracture compressibility coefficient; Effective stress; The initial effective stress before the start of extraction; The matrix shrinkage coefficient; To solve the absorption of contraction strain; Described by Langmuir-type equations; Using measured permeability data under different fracturing conditions, the model parameters were inverted and calibrated, and the accuracy was verified.

6. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, Step five specifically includes: By combining the control equations of coal and rock deformation field, gas seepage field, and gas desorption and diffusion equation, a multi-field coupled mathematical model of stress and seepage is constructed for the entire process of coal seam floor fracturing and extraction. Based on the actual geological conditions of the mine, an engineering-scale geometric model of coal seam-floor strata is established, and the fracture network of the optimized fracturing scheme is embedded into the geometric model. The ground stress boundary, gas pressure boundary, and extraction borehole negative pressure boundary are set. The multi-field coupled model is solved using the finite element method, and the spatiotemporal distribution data of coal seam stress field, gas pressure field, and permeability field during the extraction process are output. Based on the critical value of outburst risk, the outburst reduction zone is dynamically identified and divided, and the outburst reduction range and the time to achieve the outburst reduction target are output for different extraction times.

7. The method for topological quantification and permeability enhancement / outburst reduction prediction of hydraulic fracturing fractures in the coal seam floor as described in claim 1, characterized in that, Step six specifically includes: By changing individual fracturing and extraction parameters while keeping other parameters constant, the sensitivity of each parameter to the impact on the peak burst reduction period and range is analyzed by comparing simulation results, and key control parameters are identified. Parameter scanning simulations are then performed on the key control parameters within their value ranges, and the simulation results are correlated with the input parameters to create a three-dimensional quantitative relationship map of fracturing parameters, fracture topology, and peak burst reduction period. Based on the engineering requirements for the peak burst reduction range and period, the optimal fracturing process parameters and extraction system layout are selected through reverse lookup using the quantitative map, forming an application scheme for bottom plate fracturing permeability enhancement and peak burst reduction engineering technology.

Images