Fractured shale reservoir fluid injection fracture activation process simulation method

Abstract

The invention discloses a fractured shale reservoir fluid injection fracture activation process simulation method, which belongs to the technical field of unconventional oil and gas exploitation, and comprises the following steps: inputting a multi-scale fracture in-situ parameter characteristic data set into an inversion model, and performing inversion on multi-scale fracture in-situ characteristic parameters to obtain an inversion result; constructing a three-dimensional multi-scale fracture geometric model; an inversion result is fused into a geometric model, and a three-dimensional multi-scale discrete fracture network model is constructed; based on the shale macroscopic rock mechanical parameter distribution function, constructing a macroscopic heterogeneous shale reservoir matrix model; embedding the three-dimensional multi-scale discrete fracture network model into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model; according to the shale reservoir model, a fractured shale reservoir fluid injection fracture activation model is constructed, and the fracture activation slippage process is subjected to simulation characterization. The fracture activation process numerical model constructed by the method can accurately describe the dynamic evolution process of the fracture from a stable state to instability activation.

Description

Technical Field

[0001] This invention relates to the field of unconventional oil and gas extraction technology, and more specifically to a method for simulating the fluid injection fracture activation process in fractured shale reservoirs. Background Technology

[0002] With the deepening of global energy structure transformation strategies, the efficient and green development of unconventional oil and gas resources, as well as gas storage and carbon storage technologies, have become core issues in the energy and environmental protection fields. Shale reservoirs, as important unconventional oil and gas resources, rely heavily on the complex fracture networks formed by large-scale hydraulic fracturing technology for their commercial development. While injecting external fluids into fractured shale reservoirs can effectively modify and utilize reservoir oil and gas and improve conductivity, it often activates fault slip and generates severe casing deformation risks, significantly impacting single-well production.

[0003] However, current methods for characterizing fracture properties at multiple scales, particularly the activation of fractures during fracturing fluid injection, are often idealized and homogenized, failing to meet the urgent needs for technological breakthroughs and efficient development of complex deep unconventional energy resources. Against this backdrop, a technique has emerged that utilizes deep learning to construct multi-scale fractured formations and employs multi-material field methods to characterize in-situ reservoir fracture activation and slip processes.

[0004] Since the fracturing and seepage process of injecting external fluids into fractured shale reservoirs is a complex nonlinear problem involving strong coupling of multiple physical fields such as heat, fluid, solid, and damage, traditional methods for injecting fluids into fractured shale reservoirs mainly use numerical simulation methods. These methods often simplify or idealize the fracture network, fail to consider the differences in physical properties of multi-level fractures, and are difficult to truly reflect the discrete fracture system under in-situ geological conditions, resulting in significant deviations in the simulation of fracture activation dynamics. Summary of the Invention

[0005] To address the problems existing in the above-mentioned fields, this invention proposes a method and system for simulating the fracture activation process of fluid injection in fractured shale reservoirs. The constructed fracture activation model of fluid injection in fractured shale reservoirs clarifies the direct correlation between stress changes and temperature and pore pressure changes, and can accurately characterize the dynamic evolution process of fractures from a stable state to unstable activation.

[0006] To address the aforementioned technical problems, this invention discloses a method for simulating the fluid injection fracture activation process in fractured shale reservoirs, comprising the following steps: Acquire well logging imaging data, microseismic data, and seismic data of shale reservoirs in the study area, and construct a multi-scale in-situ fracture parameter feature dataset. The in-situ parameter feature dataset of multi-scale fractures is input into a pre-trained inversion model to invert the in-situ feature parameters of multi-scale fractures and obtain the inversion results. The inversion model is based on a physical information neural network and adds a residual loss term coupled with the shale seepage law and mechanical strength constraint after the output layer. Construct a three-dimensional multi-scale crack geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete crack network model; Based on the distribution function of macroscopic rock mechanical parameters of shale, a macroscopic heterogeneous shale reservoir matrix model is constructed; a three-dimensional multi-scale discrete fracture network model is embedded into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and the heterogeneity of the reservoir. Based on the shale reservoir model, a fracture activation criterion is established using Coulomb's friction law. Taking into account the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics coupling effect of heat-fluid-solid-damage caused by external fluid injection, a fracture activation model for fractured shale reservoirs is constructed to simulate and characterize the fracture activation and slip process.

[0007] Preferably, the step of inputting the multi-scale fracture in-situ parameter feature dataset into a pre-trained inversion model to invert the multi-scale fracture in-situ feature parameters and obtain the inversion result; wherein, the inversion model is based on a physical information neural network, and a residual loss term coupled with the shale seepage law and mechanical strength constraint is added after the output layer, specifically including: The input parameters of the physical information neural network are determined as follows: The output parameters are obtained as follows: ; in, The first physical information neural network represents the... One output item, The first physical information neural network represents the... One input item, i Indexes representing in-situ characteristic parameters of multi-scale cracks. Indicates the first The weighting coefficients corresponding to the in-situ characteristic parameters of multi-scale cracks Indicates the bias term; n The dimension of the dataset representing the in-situ feature parameters of multi-scale cracks; The introduced residual loss term couples the shale seepage law and mechanical strength constraints, where: The residual loss term under the seepage law constraint is: ; in, Indicates the permeability of the crack. Indicates the first Pressure gradients corresponding to in-situ characteristic parameters of multi-scale cracks; N The number of items indicates the total number of data items in the entire set. The residual loss term for the mechanical strength constraint is: ; in, and These represent the first and second digits of the shale reservoir. The in-situ characteristic parameters of multi-scale cracks, corresponding to shear stress and normal stress. and This indicates the predicted cohesion and internal friction angle of the crack; Based on the residual loss terms constrained by the seepage law and the residual loss terms constrained by the mechanical strength, the total loss function is defined as follows: ; in, Indicates hyperparameters, This represents the prediction loss of a Bayesian deep learning network; The output parameters of the physical information neural network are optimized by the total loss function to obtain the inversion results, including mechanical parameters and seepage parameters under multi-scale crack differences. The mechanical parameters are the shear stress and normal stress of the shale, the predicted cohesion and internal friction angle of the fracture; the seepage parameters are the permeability and pressure gradient of the fracture.

[0008] Preferably, the construction of a three-dimensional multi-scale crack geometric model and the fusion of the inversion results into the geometric model to construct a three-dimensional multi-scale discrete crack network model specifically includes: Based on the MATLAB platform, a three-dimensional multi-scale crack geometric model is constructed by randomly generating multiple sets of disk-shaped cracks. The constructed three-dimensional multi-scale crack geometry model is divided into three-dimensional meshes, and the inversion results are assigned to each mesh node. Each mesh point is assigned corresponding rock mechanics parameters to simulate the rock mechanics properties of the cracks, thus constructing a three-dimensional multi-scale discrete crack network model with in-situ parameter features.

[0009] Preferably, the construction of a macroscopic heterogeneous shale reservoir matrix model based on the distribution function of macroscopic rock mechanical parameters of shale specifically includes: Based on the Weibull distribution, the distribution function of macroscopic rock mechanical parameters of shale is constructed as follows: ; in, These represent rock mechanical parameters, the average value of rock mechanical parameters, and the heterogeneity coefficient, respectively. Based on the distribution function of macroscopic rock mechanical parameters of shale, a heterogeneity coefficient is introduced to simulate the spatial distribution of rock mechanical parameters of shale reservoir, thereby generating a heterogeneous shale reservoir geological model and constructing a macroscopic heterogeneous shale reservoir model.

[0010] Preferably, based on the shale reservoir model, a fracture activation criterion is established using Coulomb's friction law. This criterion comprehensively considers the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics coupling effect of heat-fluid-solid-damage caused by external fluid injection. A fracture activation model for fractured shale reservoirs is constructed to simulate and characterize the fracture activation and slip process. Specifically, this includes: Based on the macroscopically heterogeneous shale reservoir model, multiple sets of shale reservoir model parameters were obtained; Using multiple sets of shale reservoir model parameters as input parameters for the block discrete element development platform, we constructed a mechanical and seepage property parameter model considering different fracture scales, and assigned corresponding mechanical and seepage properties to fractures of different scales. The discrete element method is used to describe the deformation of fracture surfaces and reservoir rocks. The reservoir is discretized into tetrahedral meshes, and a linear elastic constitutive model is used to describe the deformation of the rock matrix. Stress field equations are then constructed. The fluid flow within the fracture is described based on the simplified Navier-Stokes equations, the flow field equations are constructed, and the relationship between fluid pressure and fracture volume change is established. Based on Fourier's law and the energy conservation equation, the heat exchange between fluids and solids is described, and the temperature field equation is constructed. Based on the stress field equation, flow field equation and temperature field equation, the thermal strain increment and pore pressure change caused by thermal expansion are introduced to construct a thermal-fluid-solid coupling mechanism, and the thermal stress increment and pore pressure increment are coupled into the mechanical constitutive equation to establish the thermal-fluid-solid-damage multiphysics field coupling control equation. Using Coulomb's law of friction as the basic model of fracture, a crack activation criterion is established. When the normal stress exceeds the tensile strength, it is judged as tensile failure; when the shear stress exceeds the maximum allowable shear force based on cohesion, internal friction angle and normal stress, it is judged as shear failure. The incremental normal and tangential stresses caused by thermal stress are superimposed on Coulomb's friction law to correct the critical condition for crack activation; once crack activation is determined, tensile strength and cohesion are set to zero to simulate the damage evolution process. Based on mechanical and seepage properties, and using the multi-physics field coupling control equations of heat-fluid-solid-damage, a fracture activation model for fluid injection in fractured shale reservoirs is constructed. Numerical simulation of the fracture activation and slip process is performed to obtain multi-scale fracture activation and propagation slip maps during fracturing injection.

[0011] Preferably, the thermo-fluid-solid-damage multiphysics coupling control equation is: ; in, and It is the first step in the fluid injection process. i The crack to the first j The stress tensor components and strain tensor components corresponding to each crack; This is the Kronecker operator notation, and the table below represents the ij-th component of a symmetric tensor matrix, where... ; b It is the Biot coefficient; The shear modulus of the rock; The bulk modulus of the rock. , E The elastic modulus of the rock. v The Poisson's ratio of the rock; It is the linear thermal expansion coefficient of the rock; This is the sum of the strain tensors of the rock in different directions.

[0012] Preferably, the simulation and characterization of the fracture activation and slip process further includes determining the fracture activation distribution characteristics based on the obtained multi-scale fracture activation and propagation slip maps during hydraulic fracturing injection, combined with microseismic event evaluation, specifically including: According to the clustering mechanism, when two slip segments in the multi-scale crack activation and extension slip map are adjacent to each other and the slip duration overlaps, the slip nodes within the preset radius R are merged and the microseismic slip segments are updated. The slip velocity is determined based on the updated microseismic slip segment. When the slip velocity is greater than the velocity threshold, it is determined to be a seismic slip; otherwise, it is determined to be a quiescent slip. For the slip segment determined to be seismically slipped, the seismic moment and magnitude are calculated based on the slip amount of the crack, the slip fracture area and the rock shear modulus, and microseismic spheres are generated. Based on the spatial location of microseismic spheres, a one-to-one correspondence between fault slip rate and microseismic events is established to determine the spatial distribution characteristics of fault activation.

[0013] Preferably, the acquisition of well logging imaging data, microseismic data, and seismic data of shale reservoirs in the study area, and the construction of a multi-scale in-situ fracture parameter feature dataset, specifically includes: Collect well logging imaging data, microseismic data, and seismic data from oil and gas field construction sites; Grayscale and binary processing were performed on the logging imaging data to obtain in-situ characteristic parameters of well perimeter microfractures at the micrometer to centimeter level; Based on microseismic data and the in-situ characteristic parameters of meter-scale cracks were obtained using the inversion solution derived from the focal mechanism. The in-situ characteristic parameters of 100-meter-level fractures were obtained by processing seismic data using the ant body method. By fusing in-situ feature parameters of micro-fractures, micro-fractures, and fractures at the micrometer to centimeter scale, as well as those at the meter and hundred-meter scale, a multi-scale fracture in-situ parameter feature dataset is constructed.

[0014] Preferably, it also includes a simulation system for the fluid injection fracture activation process in fractured shale reservoirs, comprising: A multi-scale dataset construction module is used to acquire well logging imaging data, microseismic data and seismic data of shale reservoirs in the study area, and to construct a multi-scale fracture in-situ parameter feature dataset. The parameter inversion module is used to input the multi-scale fracture in-situ parameter feature dataset into the pre-trained inversion model, invert the multi-scale fracture in-situ feature parameters, and obtain the inversion results; wherein, the inversion model is based on a physical information neural network, and a residual loss term coupled with shale seepage law and mechanical strength constraint is added after the output layer; The shale reservoir model construction module is used to construct a three-dimensional multi-scale fracture geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete fracture network model; construct a macroscopic heterogeneous shale reservoir matrix model based on the distribution function of macroscopic rock mechanical parameters of shale; and embed the three-dimensional multi-scale discrete fracture network model into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and reservoir heterogeneity. The fracture activation and slip process simulation module is used to establish a fracture activation criterion based on the shale reservoir model and the Coulomb friction law. It also comprehensively considers the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics field coupling effect of heat-fluid-solid-damage caused by external fluid injection, and constructs a fracture activation model of fractured shale reservoirs by fluid injection to simulate and characterize the fracture activation and slip process.

[0015] Compared with the prior art, the present invention has the following beneficial effects: The proposed method for simulating the fluid injection and fracture activation process in fractured shale reservoirs utilizes an inversion model that combines data-driven and physical constraints. This approach ensures that the inversion results better reflect physical laws, significantly improving the accuracy of in-situ inversion of multi-scale fracture characteristic parameters compared to purely data-driven models. This provides fundamental parameter support for subsequent fracture activation simulations that closely align with actual reservoir conditions. The shale reservoir model accurately reflects the multi-scale discrete fracture network system and reservoir heterogeneity, making the model more closely resemble in-situ geological conditions. The constructed fracture activation model for fractured shale reservoirs by fluid injection comprehensively considers heat conduction, thermal expansion, fluid seepage, solid deformation, multi-field coupling, and damage evolution. This synergistic mechanism of heat-fluid-solid-damage multi-physics fields overcomes the shortcomings of existing numerical simulation methods that simplify or idealize fracture networks. It can accurately characterize the fracture instability process, accurately predict the critical conditions for fracture activation under fluid injection, and for the first time clarify the direct correlation between stress change and temperature and pore pressure changes. It makes up for the shortcomings of the traditional Coulomb criterion, which only considers mechanical factors and ignores the thermal stress increment and pore pressure change caused by fluid injection. It can finely characterize the coupling mechanism of temperature field-seepage field-stress field during fluid injection, accurately depict the dynamic evolution process of fractures from a stable state to instability activation, and solve the technical pain point of uncontrolled fault activation caused by blind injection in unconventional oil and gas development. Attached Figure Description

[0016] Figure 1 This is a flowchart of the simulation method for the fluid injection fracture activation process in fractured shale reservoirs proposed in this invention; Figure 2 This is a well logging imaging data diagram provided in an embodiment of the present invention; Figure 3 A flowchart of the source mechanism inversion solution provided in an embodiment of the present invention; Figure 4 This is a prediction diagram of cracks in the ant body provided in an embodiment of the present invention; Figure 5 A structural diagram of a Bayesian deep learning framework for seepage theory and mechanical strength constraints provided in an embodiment of the present invention; Figure 6 This is a diagram showing the elastic modulus distribution of a shale reservoir model provided in an embodiment of the present invention. Figure 7 A diagram of a three-dimensional multi-scale discrete crack network model provided in an embodiment of the present invention; Figure 8 Numerical model for fluid injection in multi-scale fractured shale reservoirs provided for embodiments of the invention; Figure 9 The figure shows the simulation results of the pore pressure increment of fluid injection in multi-scale fractured shale reservoirs provided in the embodiments of the present invention; Figure 10This is a simulation result of the fracture opening displacement of fluid injection in a multi-scale fractured shale reservoir provided in an embodiment of the present invention; Figure 11 The figure shows the simulation results of the fracture slip rate of fluid injection in multi-scale fractured shale reservoirs provided in the embodiments of the present invention; Figure 12 The figure shows the simulation results of fluid injection fracture slip in multi-scale fractured shale reservoirs provided in the embodiments of the present invention. Figure 13 This is a spatiotemporal combination characterization mechanism diagram of adjacent sliding nodes in a microseismic event provided in an embodiment of the present invention; Figure 14 The figure shows the simulation results of multi-scale fractured shale reservoir fluid injection fracture-activated microseismic event points provided in the embodiments of the present invention. Detailed Implementation

[0017] The following will refer to the appendices in the embodiments of the present invention. Figures 1-14 The technical solutions in the embodiments of the present invention will be clearly and completely described. It should be understood that the terminology used in the present invention is only for describing particular implementation methods and is not intended to limit the present invention.

[0018] Example This embodiment selects well L2 in block xx of the shale reservoir as the study area, such as Figure 1 As shown, this invention proposes a method for simulating the fluid injection fracture activation process in fractured shale reservoirs, comprising the following steps: S1: Collect logging imaging data, microseismic data, and seismic data from well L2 (the construction site of the oil and gas field) in block xx.

[0019] Grayscale and binary processing were performed on well logging imaging data to obtain in-situ characteristic parameters of wellbore microfractures at the micrometer to centimeter scale, such as... Figure 2 The image shown is a diagram of in-situ characteristic parameters of microfractures obtained from well logging imaging data.

[0020] Based on microseismic data and focal mechanism inversion solutions, in-situ characteristic parameters of meter-scale cracks were obtained, such as... Figure 3 The diagram shows the flow chart of the focal mechanism inversion solution. Microseismic event points are mainly acquired by the detector. The spatial location points are obtained through feature extraction and mapping. The focal mechanism inversion solution can be obtained by using focal mechanism inversion. The crack features are characterized by obtaining the beach ball. Finally, the b-value feature is used to screen natural cracks, and the fracture surface is extracted by weighted three-dimensional plane fitting, forming a crack model based on microseismic events.

[0021] The ant-body method is used to process seismic data and obtain in-situ characteristic parameters of faults at the hundred-meter level, such as... Figure 4 The image shown is a prediction diagram of cracks in the ant's body.

[0022] Finally, the in-situ feature parameters of micro-fractures, meter-scale fractures, and hundred-meter-scale fractures around the well are fused to construct a multi-scale fracture in-situ parameter feature dataset.

[0023] In step S1, this method integrates well logging imaging (micrometer-centimeter level), microseismic (meter level), and seismic (hundred-meter level) data to construct a multi-scale fracture parameter set covering six orders of magnitude. This full-scale characterization avoids the limitations of single-scale models. For example, traditional methods may ignore the dominant role of microfractures in seepage or the macroscopic influence of fracture zones on stress distribution.

[0024] Outlier removal and data normalization are performed on the constructed dataset to obtain a preprocessed dataset; then, the dataset is divided into training and test sets according to a preset ratio.

[0025] S2: Input the dataset into the inversion model to obtain the mechanical parameters and seepage parameters under the multi-scale fracture differences through inversion; wherein, the inversion model is based on the physical information neural network, and a residual loss term coupled with the shale seepage law and mechanical strength constraint is added after the output layer.

[0026] like Figure 5 As shown, the physical information neural network provided in this embodiment adopts a Bayesian deep learning network. The original Bayesian deep learning network includes an input layer, a hidden layer, and an output layer. The inversion model constructed in this invention is based on the Bayesian deep learning network, with physical constraints added after the output layer. The physical constraints include a residual loss term that couples the shale seepage law and mechanical strength constraints.

[0027] According to Bayes' theorem, the input parameters of a Bayesian deep learning network are determined as follows: The output parameters are obtained as follows: ; in, The first term of a Bayesian deep learning network is represented by the second term. One output item, The first term of a Bayesian deep learning network is represented by the second term. One input item, i Indexes representing in-situ characteristic parameters of multi-scale cracks. Indicates the first The weighting coefficients corresponding to the in-situ characteristic parameters of multi-scale cracks Indicates the bias term; n The dimension of the dataset representing the in-situ feature parameters of multi-scale cracks; The constructed inversion model incorporates physical constraints coupled with residual loss terms from shale seepage laws and mechanical strength constraints, where: The residual loss term under the seepage law constraint is: ; in, Indicates the permeability of the crack. Indicates the first Pressure gradients corresponding to in-situ characteristic parameters of multi-scale cracks; N The number of items indicates the total number of data items in the entire set. The residual loss term for the mechanical strength constraint is: ; in, and These represent the first and second digits of the shale reservoir. The in-situ characteristic parameters of multi-scale cracks, corresponding to shear stress and normal stress. and This indicates the predicted cohesion and internal friction angle of the crack; Based on the residual loss terms constrained by the seepage law and the residual loss terms constrained by the mechanical strength, the total loss function is defined as follows: ; in, Indicates hyperparameters, This represents the prediction loss of a Bayesian deep learning network.

[0028] By optimizing the output parameters of the Bayesian deep learning network using the total loss function, mechanical and seepage parameters under multi-scale fracture differences are obtained. Among them, the mechanical parameters are the shear stress and normal stress of shale, the cohesion and internal friction angle predicted by the fracture, and the seepage parameters are the permeability and pressure gradient of the fracture.

[0029] The inversion model is trained using the training set to obtain the trained inversion model, and its performance is evaluated on the test set.

[0030] The trained inversion model can stably invert the in-situ characteristic parameters of multi-scale cracks. Based on the trained inversion model, the in-situ characteristic parameters of multi-scale cracks are inverted to obtain the inversion results.

[0031] The inversion model constructed in this invention introduces a residual loss term that couples seepage theory with mechanical strength based on a Bayesian deep learning network. By combining data-driven and physical constraints, it significantly improves the accuracy of the inversion results and provides key parameter inputs for fracture activation simulation.

[0032] The mean absolute percentage error (MAPE) of the inversion model on the training set and the test set were 8.76% and 13.9%, respectively, which verifies that the inversion model has high accuracy in the inversion of in-situ characteristic parameters of multi-scale cracks.

[0033] S3: Construct a three-dimensional multi-scale crack geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete crack network model.

[0034] Based on the MATLAB platform, a three-dimensional multi-scale fracture geometric model is constructed by randomly generating multiple sets of disk-shaped fractures. The constructed three-dimensional multi-scale fracture geometric model is divided into three-dimensional meshes, and the mechanical parameters and seepage parameters obtained by inversion are assigned to each mesh node. Each mesh point is assigned corresponding rock mechanical parameters to simulate the rock mechanical properties of the fracture, thereby constructing a three-dimensional multi-scale discrete fracture network model with in-situ parameter features and extracting multi-scale fracture network features.

[0035] In step S4, based on the Weibull distribution, the distribution function of macroscopic rock mechanical parameters of shale is constructed as follows: ; in, These represent rock mechanical parameters, the average value of rock mechanical parameters, and the heterogeneity coefficient, respectively. S4: Based on the distribution function of macroscopic rock mechanical parameters of shale, a macroscopic heterogeneous shale reservoir matrix model is constructed; a three-dimensional multi-scale discrete fracture network model is embedded into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and the heterogeneity of the reservoir.

[0036] Based on the distribution function of macroscopic rock mechanical parameters of shale, a heterogeneity coefficient is introduced to simulate the spatial distribution of rock mechanical parameters of shale reservoirs, generating a heterogeneous shale reservoir geological model. This allows for the construction of a shale reservoir model that can realistically reflect the multi-scale discrete fracture network system and reservoir heterogeneity, serving as a macroscopic heterogeneous shale reservoir model to obtain multiple sets of shale reservoir models.

[0037] Specifically, in this embodiment, a three-dimensional multi-scale discrete fracture network model containing in-situ parameter features is embedded into a distribution function to construct a shale reservoir geological model with a heterogeneity coefficient m=4.8, forming a shale reservoir model that can realistically reflect the multi-scale discrete fracture network system and reservoir heterogeneity; such as Figure 6 The figure shows the distribution of elastic modulus in a shale reservoir model. The spatial distribution of shale reservoir rock mechanical parameters is simulated using the distribution function of macroscopic rock mechanical parameters. Embedding a discrete fracture network into the macroscopic heterogeneous shale reservoir matrix model can accurately reflect the true characteristics of shale reservoirs.

[0038] The distribution function of macroscopic rock mechanical parameters of shale constructed by the Weibull distribution can simulate the spatial variation of parameters such as reservoir permeability and porosity.

[0039] Embedding discrete fracture networks into the distribution function of macroscopic rock mechanical parameters of shale allows for the analysis of local response differences in fracture activation. For example, fractures propagate faster in high-permeability zones and require higher injection pressure in low-permeability zones, thereby guiding the optimization of fracturing segment spacing and displacement.

[0040] S5: Based on the shale reservoir model, a fracture activation criterion is established using Coulomb's friction law. Taking into account the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics field coupling effect of heat-fluid-solid-damage caused by external fluid injection, a fracture activation model for fractured shale reservoirs is constructed to simulate and characterize the fracture activation and slip process.

[0041] S51: The discrete element method is used to describe the deformation of fracture surfaces and reservoir rocks. The reservoir is discretized into tetrahedral meshes, and a linear elastic constitutive model is used to describe the deformation of the rock matrix, thus constructing the stress field equations.

[0042] Specifically, this invention uses the fixed points of each tetrahedral mesh as computational nodes, and determines the motion equations for each node as follows: ; in, It is the surface area of ​​the rock. It is the unit normal vector of the surface. It is the resultant force of all external forces corresponding to the i-th crack. It is the gravitational acceleration corresponding to the i-th crack. Let be the displacement corresponding to the i-th crack.

[0043] Specific F i The sum of the three forces at the i-th crack point on the grid is expressed as: ; in, It is the external load of the i-th crack. This represents the contact force of the i-th crack, which exists only at mesh points along the block boundary. The internal stress of the i-th crack near the node... The contribution is calculated as follows:

[0044] ; in, It is the stress tensor component corresponding to the i-th to j-th cracks in the fluid injection process; is the unit normal vector of the surface corresponding to the j-th crack; z is the superscript of the internal stress near the node.

[0045] The sum of the total nodal stress tensors on each mesh needs to be calculated, i.e., the sum of the stress tensors corresponding to the i-th fracture to the j-th fracture during the fluid injection process is ≥ This tensor incorporates the effects of external loads and gravity. Gravity The calculation method is as follows:

[0046] ; in, m This represents the quality of fractured shale; This is the acceleration due to gravity.

[0047] When the object is in equilibrium, the sum of the total nodal stress tensors on each mesh. >0; otherwise, according to the finite difference form of Newton's second law of motion, the total acceleration of each node on each grid can be obtained as follows: ; The superscript indicates the time corresponding to the calculated variable. In each time step, the strain tensor between any two nodes on each grid is... Rotation angle tensor The relationship between the displacement tensor and the displacement tensor can be expressed as:

[0048] ; in, Represents nodes during fluid injection. i To the node j The displacement tensor, Represents nodes during fluid injection. j To the node i The displacement tensor.

[0049] like Figure 7 The figure shown is a three-dimensional multi-scale discrete crack network model with in-situ parameter features constructed in this invention, which can accurately describe the geometric shape (such as length, width, and aperture) and spatial distribution (such as connectivity and azimuth) of cracks.

[0050] This invention employs a linear elastic constitutive model to describe the deformation of the rock matrix, constructs the stress field equations, and determines the form of the stress increment as follows: ; in, and It is Lamé's constant. This represents the increment of the stress tensor components corresponding to the i-th fracture to the j-th fracture during the fluid injection process in fractured shale. This represents the increment of the strain tensor components corresponding to the i-th fracture to the j-th fracture during the fluid injection process in fractured shale. It is the increment of volumetric strain. It is the Kronecker function.

[0051] S52: The fluid flow within the fracture is described based on the simplified Navier-Stokes equations. The flow field equations are constructed, and the relationship between fluid pressure and fracture volume change is established.

[0052] This invention governs the fluid flow in a crack according to the simplified Navier-Stokes equations. When the fluid flows between two nearly parallel impermeable boundaries, the flow field equations are simplified to the Reynolds equations:

[0053] ; in, It is a point on a plane The distance between impermeable boundaries; The water head corresponding to the i-th crack; Represents gravitational acceleration; For fluid density, For fluid viscosity, For height; It is the pressure of the fluid, the integral equation This represents the distance between impermeable boundaries. For fracture surfaces between rock blocks, the integral range is the distance between two fracture surfaces, and the fluid velocity corresponding to the fracture width of the i-th fracture is:

[0054] ; in, It is the permeability of a single crack. is the hydraulic conductivity coefficient.

[0055] Fluid flow in rock reservoirs is a hydraulic-mechanical coupling problem; fluid pressure within fractures affects normal fracture deformation. Walsh derived that the pore space volume does not change with normal stress. and pore pressure The expression for the effective stress under varying conditions is:

[0056] ; in, It is the volume of the pore space. It is the compressibility of the rock surrounding the crack. In addition, effective stress... The equation can be expressed as:

[0057] ; The equation can be simplified to: ; Pore ​​pressure also affects the strength of rock fractures through the effective stress principle. When using effective stress to assess the tensile or shear strength of rock fractures, the width of a typical hydraulic fracture can be expressed as:

[0058] ; in, Indicates the crack aperture under zero normal stress. This represents the increment of normal displacement of the joint / fracture. Considering the net flow rate into the fracture and the volumetric mechanics caused by the movement of the surrounding rock mass, it is necessary to update the fluid pressure within the fracture. The new fluid pressure is:

[0059] ; in: ; in, It is the fluid pressure at the previous time step. It is the bulk modulus of the fluid, here and These are the crack volumes at the current and previous time steps, respectively. and These are the inflow rates of the i-th crack, the adjacent crack, and the source term, respectively. and These are the fluid densities under injection and formation conditions, respectively.

[0060] S53: Based on Fourier's law and the energy conservation equation, describe the heat exchange between fluid and solid, and construct the temperature field equation.

[0061] For the thermal convection problem during flow within a crack, heat transfer between fluids within the crack follows Fourier's law, while heat transfer between the fluid and the solid is a forced convection heat transfer problem. The energy conservation equation for the fluid is: ; ; in, It is fluid density. It is the specific heat of the fluid. It is the specific heat flux in a fluid. It indicates the change in fluid velocity. It indicates the amount of temperature change of the fluid. The contact area per unit fluid. The heat transfer coefficient between the fluid and the rock. and The temperatures of the fluid and the solid block are respectively. denoted as , where is the thermal conductivity of the fluid.

[0062] S54: Based on the stress field equation, flow field equation, and temperature field equation, the thermal strain increment and pore pressure change caused by thermal expansion are introduced to construct a thermal-fluid-solid coupling mechanism, and the thermal stress increment and pore pressure increment are coupled into the mechanical constitutive equation to establish the thermal-fluid-solid-damage multi-physics field coupling control equation.

[0063] Fluid-solid-thermal coupling: Fluid pressure affects rock mass deformation, and changes in joint aperture react upon fluid flow; thermal expansion (thermal stress) caused by temperature changes leads to changes in joint aperture; heat exchange between fracturing fluid and reservoir, and heat exchange generated by friction; temperature-dependent fluid viscosity. Initial pore pressure from shale gas (methane) is also considered.

[0064] ①Thermal strain-mechanical coupling For thermo-mechanical coupling processes, thermal volumetric strain is introduced into incremental mechanics and fluid constitutive laws to account for the coupling of thermal stress and thermal pore pressure. Solving the thermal stress problem requires re-establishing the incremental stress-strain relationship, which is achieved by subtracting the portion caused by temperature change from the total strain increment.

[0065] Since free thermal expansion does not cause angular deformation in isotropic materials, the shear strain increment is unaffected and correlates with the temperature increment. The corresponding thermal strain increment related to free expansion has the following form: ; in, It is the linear thermal expansion coefficient of solid rock. It is the change in temperature. It is the Kronecker operator.

[0066] ② Fluid-coupled thermally induced pore pressure In the absence of volume deformation, the pore pressure caused by thermal expansion can be calculated as follows: ; in, It is the bulk modulus of the fluid. It is the bulk modulus of a solid. For unit reservoir porosity, then 1- This indicates the percentage of the rock skeleton within a unit reservoir. is the linear thermal expansion coefficient of the fluid; This refers to the change in pore pressure caused by thermal expansion.

[0067] This equation simply shows that the thermal expansion of water and solids should be easily compensated for by the mechanical compression of the components. If we assume the solid phase is incompressible, then:

[0068] ; ③ Thermodynamic-pore pressure coupling For the temperature-hydraulic coupling of the mechanical constitutive equation, this invention achieves this through changes in the effective normal strain rate (i.e., the total strain rate minus the thermal normal strain rate component) and pore pressure. These changes affect the effective stress; therefore, the governing equations for the thermo-fluid-solid-damage multiphysics coupling are established as follows: ; in, and These are the stress tensor components and strain tensor components corresponding to the i-th to j-th cracks during the fluid injection process; It is the Kronecker operator symbol or Kronecker function; It is the increment of volumetric strain. b It is the Biot coefficient; The shear modulus of the rock; The bulk modulus of the rock. , E The elastic modulus of the rock. v The Poisson's ratio of the rock; This is the sum of the rock strain tensors in different directions. .

[0069] S55: Using Coulomb's law of friction as the fundamental model of fracture, establish a crack activation criterion. This invention uses Coulomb's law of friction as the fundamental model for fracture, considering both shear and tensile failure as well as shear expansion. For an unactivated fracture, the tensile normal stress limit, i.e., the maximum tensile strength, is:

[0070] ; in, It is the tensile strength at break. It is the maximum tensile strength.

[0071] Maximum permissible shear stress in the normal direction for: ; in, and These are the cohesion angle and the internal friction angle at fracture, respectively. It is the normal stress at the fracture site. Once activated, both tensile strength and cohesion become zero.

[0072] ; To more accurately reflect the in-situ reservoir environment, in addition to fluid-structure interaction, the effect of temperature must also be considered in fracture activation. Therefore, the normal stress increment is given by the change in normal displacement, and the shear stress increment is controlled by the shear displacement.

[0073] ; in, The coefficient of thermal expansion at fracture; and These are the normal and tangential moduli of fracture, respectively. and These are the normal and tangential feature lengths (taken from the element size). and The normal stress increment and the shear stress increment are represented respectively. and These represent the opening displacement and the shear displacement, respectively.

[0074] S56: Based on mechanical and seepage properties, and using the multi-physics field coupling control equation of heat-fluid-solid-damage, a fracture activation model for fluid injection in fractured shale reservoirs is constructed. Numerical simulation of the fracture activation and slip process is performed to obtain multi-scale fracture activation and propagation slip maps during fracturing injection.

[0075] The numerical model of the fracture activation process of fluid injection in fractured shale reservoirs constructed in this invention is as follows: Figure 8 As shown, the model not only incorporates the dual mechanical and seepage properties of multi-scale fractures but also the heterogeneous mechanical characteristics of the shale matrix. Based on this, a fracture activation model for fluid injection in fractured shale reservoirs is constructed. For this model, the theoretical iterative calculation comprehensively considers heat conduction (e.g., heat exchange between low-temperature external fluids and high-temperature rock mass), fluid seepage (e.g., the high diffusivity of the injected fluid), thermal stress (e.g., changes in thermal stress in the rock matrix and Coulomb stress in the fracture caused by low-temperature external fluids), solid deformation (e.g., changes in the geostress field), and damage evolution (e.g., abrupt changes in permeability due to fracture propagation). This multi-field synergy can accurately predict the critical conditions for fracture activation. For example, the greater the geostress difference, the easier it is for fracture propagation to deflect along the direction of the maximum principal stress, forming a complex fracture network; the higher the fracture permeability, the easier it is for fracture extension to deflect along the fracture direction, forming a cross-cutting fracture network.

[0076] like Figure 9-12 The figure shows the simulation results of thermal-fluid-solid damage in shale reservoirs after external fluid injection, using a fracture activation model based on fluid injection into fractured shale reservoirs, as described in this invention. Figure 9 The calculation results represent the increase in pore pressure caused by simulated fluid injection into the formation and subsequent fracture, characterizing the distribution of fluid inflow into large faults and small fractures. Figure 10 The result represents the calculated fracture opening deformation in the simulation. Its magnitude shows the real-time fracture width after fracture activation, which will have a significant impact on the flow velocity of the fluid. Figure 11The result represents the calculated fracture shear rate, which reflects the degree of deformation after fracture activation. This result can be used for microseismic event synthesis and computational mapping. Figure 12 The result represents the calculated fracture shear deformation. Its magnitude shows the amount of dislocation on both sides of the fracture after the fracture is activated. If the horizontal well casing passes through the fault, the casing deformation can be calculated based on the fracture shear slip.

[0077] S57: Based on the multi-scale fracture activation and propagation slip maps obtained during the fracturing injection process, the characteristics of fracture activation distribution are determined in conjunction with microseismic event evaluation.

[0078] Microseismic event synthesis By analyzing the interaction at the contact points, the amount of crack slip can be determined to estimate the synthetic seismic event. Seismic moment ( ) and magnitude ( Calculate according to the following formula:

[0079] ; in, The shear modulus of the rock. The area of ​​each slip fracture, This indicates the amount of slip for each fracture.

[0080] Microseismic event points are mainly divided into two categories: identification and merging of identical event points, and discrimination between seismic and quiescent slip.

[0081] Based on the clustering mechanism, when two slip segments in the multi-scale crack activation and extension slip map are adjacent to each other and their slip durations overlap, the slip nodes within a preset radius R are merged, and the slip amount is merged with the preset radius range to update the microseismic slip segment. This process then allows for further identification, merging, and iterative updating of the microseismic slip segment. The specified radius R is set to 1.5 times the average side length of a tetrahedron.

[0082] like Figure 13The diagram illustrates the spatiotemporal combination mechanism of adjacent slip nodes in a microseismic event. The left diagram shows the merging process of grid slip nodes activated by the fault plane during their duration, while the right diagram shows the spatial merging process of grid slip nodes activated by the fault plane. A comprehensive analysis of these two diagrams yields the theoretical process of spatiotemporal merging of microseismic event points. In the left diagram, if node 2 and node 1 have overlapping durations during their slip time, these two nodes can merge. Similarly, if node 3 and node 1 have overlapping durations, they can also merge, ultimately forming a large-scale microseismic event that merges nodes 1, 2, and 3. The slip rate of this large-scale microseismic event during this time period can be calculated using a numerical model system. Combined with the seismic / seismic slip discrimination mode, the large-scale microseismic events with seismic activity are characterized using microseismic spheres. The slip discrimination of microseismic events primarily relies on the intensity of the slip. Based on the merged and updated microseismic slip segment results, seismic / anecdotic slip is determined according to the slip velocity of the slip segment. If the slip velocity is >0.5mm / s, it is considered seismic slip; otherwise, it is anecdotic slip.

[0083] Based on the microseismic synthetic event discrimination, a one-to-one correspondence between the fault slip rate under fluid injection and the spatial location of microseismic events was finally obtained, such as... Figure 14 As shown in the figure, the points represent microseismic spheres at the time of composite slip of the fractured fault. The spatial location of these microseismic spheres is identified and synthesized based on the grid nodes of activated slip on the fault plane, thus reflecting the spatial distribution characteristics of fault-activated slip. According to Figure 13 The theoretical process involves continuous fluid injection during fracturing, which activates and expands the grid nodes on the fault plane. Microseismic event points continuously merge, and their size reflects the magnitude of the merged slip nodes, thus indicating the relative magnitude of the activated slip to some extent. Based on this result, the relative magnitude and spatial distribution characteristics of the activated slip can be intuitively observed.

[0084] This invention also proposes a simulation system for the fluid injection fracture activation process in fractured shale reservoirs, comprising: A multi-scale dataset construction module is used to acquire well logging imaging data, microseismic data and seismic data of shale reservoirs in the study area, and to construct a multi-scale fracture in-situ parameter feature dataset. The parameter inversion module is used to input the multi-scale fracture in-situ parameter feature dataset into the pre-trained inversion model, invert the multi-scale fracture in-situ feature parameters, and obtain the inversion results; wherein, the inversion model is based on a physical information neural network, and a residual loss term coupled with shale seepage law and mechanical strength constraint is added after the output layer; The shale reservoir model construction module is used to construct a three-dimensional multi-scale fracture geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete fracture network model; construct a macroscopic heterogeneous shale reservoir matrix model based on the distribution function of macroscopic rock mechanical parameters of shale; and embed the three-dimensional multi-scale discrete fracture network model into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and reservoir heterogeneity. The fracture activation and slip process simulation module is used to establish a fracture activation criterion based on the shale reservoir model and the Coulomb friction law. It also comprehensively considers the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics field coupling effect of heat-fluid-solid-damage caused by external fluid injection, and constructs a fracture activation model of fractured shale reservoirs by fluid injection to simulate and characterize the fracture activation and slip process.

[0085] The simulation method for the fluid injection fracture activation process in fractured shale reservoirs proposed in this invention constructs a multi-scale in-situ fracture parameter feature dataset, which can avoid the limitations of a single-scale model.

[0086] The constructed inversion model incorporates a physical information neural network that couples the seepage law and mechanical strength constraints. By combining data-driven approaches with physical constraints, the accuracy of the inversion results is significantly improved.

[0087] The multi-physics coupled control equations based on Coulomb's friction law, involving heat, fluid, solid, and damage, have for the first time clarified the direct correlation between normal stress, shear stress, temperature, and pore pressure. This overcomes the shortcomings of the traditional Coulomb criterion, which only considers mechanical factors and ignores the increase in thermal stress and changes in pore pressure caused by fluid injection. It can finely characterize the coupling mechanism of "temperature field-seepage field-stress field" during fluid injection and accurately depict the dynamic evolution of cracks from a stable state to instability activation.

[0088] A numerical model of the fracture activation process during fluid injection in fractured shale reservoirs has been developed. This model enables dynamic simulation and visualization of the timing of fault instability, slip direction, and propagation path during fluid injection, addressing the technical challenge of "uncontrolled fault activation due to blind injection" in unconventional oil and gas development. This predictive capability can identify potential large-scale fracture slip risks in advance, providing quantitative basis for optimizing fluid injection parameters (injection pressure, injection rate, and injection temperature) and effectively reducing the geological hazard risks during reservoir fracturing.

[0089] The constructed fracture activation model for fractured shale reservoirs by fluid injection comprehensively considers the thermal stress on the shale reservoir caused by external fluids, the changes in pore pressure field caused by thermal expansion, the increase in normal and tangential stress caused by thermal stress, the in-situ shale reservoir environment, and the coupling effect of multiple physical fields of heat, fluid, solid, and damage. This in-situ shale reservoir environment and the coupling effect of multiple physical fields of heat, fluid, solid, and damage can accurately predict the critical conditions for fracture activation, providing an efficient, safe, and economical solution for unconventional oil and gas development.

[0090] The above description is only a preferred embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and inventive concept of the present invention, should be covered within the scope of protection of the present invention.

[0091] Furthermore, unless otherwise stated, all technical and scientific terms used in this invention have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains. All references to this specification are incorporated by way of citation to disclose and describe methods relating to those references. In the event of any conflict with any incorporated reference, the content of this specification shall prevail.

Claims

1. A method for simulating the fracture activation process of fluid injection in fractured shale reservoirs, characterized in that, Includes the following steps: Acquire well logging imaging data, microseismic data, and seismic data of shale reservoirs in the study area, and construct a multi-scale in-situ fracture parameter feature dataset. The in-situ parameter feature dataset of multi-scale fractures is input into a pre-trained inversion model to invert the in-situ feature parameters of multi-scale fractures and obtain the inversion results. The inversion model is based on a physical information neural network and adds a residual loss term coupled with the shale seepage law and mechanical strength constraint after the output layer. Construct a three-dimensional multi-scale crack geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete crack network model; Based on the distribution function of macroscopic rock mechanical parameters of shale, a macroscopic heterogeneous shale reservoir matrix model is constructed; a three-dimensional multi-scale discrete fracture network model is embedded into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and the heterogeneity of the reservoir. Based on the shale reservoir model, a fracture activation criterion is established using Coulomb's friction law. Taking into account the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics coupling effect of heat-fluid-solid-damage caused by external fluid injection, a fracture activation model for fractured shale reservoirs is constructed to simulate and characterize the fracture activation and slip process.

2. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 1, characterized in that, The process involves inputting a multi-scale fracture in-situ parameter feature dataset into a pre-trained inversion model to invert the multi-scale fracture in-situ feature parameters and obtain the inversion results. The inversion model is based on a physical information neural network, and a residual loss term coupling the shale seepage law and mechanical strength constraints is added after the output layer. Specifically, this includes: The input parameters of the physical information neural network are determined as follows: The output parameters are obtained as follows: ； in, The first physical information neural network represents the... One output item, The first physical information neural network represents the... One input item, i Indexes representing in-situ characteristic parameters of multi-scale cracks. Indicates the first The weighting coefficients corresponding to the in-situ characteristic parameters of multi-scale cracks Indicates the bias term; n The dimension of the dataset representing the in-situ feature parameters of multi-scale cracks; The introduced residual loss term couples the shale seepage law and mechanical strength constraints, where: The residual loss term under the seepage law constraint is: ； in, Indicates the permeability of the crack. Indicates the first Pressure gradients corresponding to in-situ characteristic parameters of multi-scale cracks; N The number of items indicates the total number of data items in the entire set. The residual loss term for the mechanical strength constraint is: ； in, and These represent the first and second digits of the shale reservoir. The in-situ characteristic parameters of multi-scale cracks, corresponding to shear stress and normal stress. and This indicates the predicted cohesion and internal friction angle of the crack; Based on the residual loss terms constrained by the seepage law and the residual loss terms constrained by the mechanical strength, the total loss function is defined as follows: ； in, Indicates hyperparameters, This represents the prediction loss of a Bayesian deep learning network; The output parameters of the physical information neural network are optimized by the total loss function to obtain the inversion results, including mechanical parameters and seepage parameters under multi-scale crack differences. The mechanical parameters are the shear stress and normal stress of the shale, the predicted cohesion of the fracture, and the internal friction angle; the seepage parameters are the permeability of the fracture and the pressure gradient.

3. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 1, characterized in that, The construction of a three-dimensional multi-scale crack geometric model and the integration of inversion results into the geometric model to construct a three-dimensional multi-scale discrete crack network model specifically include: Based on the MATLAB platform, a three-dimensional multi-scale crack geometric model is constructed by randomly generating multiple sets of disk-shaped cracks. The constructed three-dimensional multi-scale crack geometry model is divided into three-dimensional meshes, and the inversion results are assigned to each mesh node. Each mesh point is assigned corresponding rock mechanics parameters to simulate the rock mechanics properties of the cracks, thus constructing a three-dimensional multi-scale discrete crack network model with in-situ parameter features.

4. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 1, characterized in that, The construction of a macroscopic heterogeneous shale reservoir matrix model based on the distribution function of macroscopic rock mechanical parameters of shale specifically includes: Based on the Weibull distribution, the distribution function of macroscopic rock mechanical parameters of shale is constructed as follows: ； in, These represent rock mechanical parameters, the average value of rock mechanical parameters, and the heterogeneity coefficient, respectively. Based on the distribution function of macroscopic rock mechanical parameters of shale, a heterogeneity coefficient is introduced to simulate the spatial distribution of rock mechanical parameters of shale reservoir, thereby generating a heterogeneous shale reservoir geological model and constructing a macroscopic heterogeneous shale reservoir model.

5. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 1, characterized in that, Based on the shale reservoir model, a fracture activation criterion is established using Coulomb's friction law. Taking into account the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics coupling effect of heat-fluid-solid-damage caused by external fluid injection, a fracture activation model for fractured shale reservoirs is constructed. The model simulates and characterizes the fracture activation and slip process, specifically including: Based on the macroscopically heterogeneous shale reservoir model, multiple sets of shale reservoir model parameters were obtained; Using multiple sets of shale reservoir model parameters as input parameters for the block discrete element development platform, we constructed a mechanical and seepage property parameter model considering different fracture scales, and assigned corresponding mechanical and seepage properties to fractures of different scales. The discrete element method is used to describe the deformation of fracture surfaces and reservoir rocks. The reservoir is discretized into tetrahedral meshes, and a linear elastic constitutive model is used to describe the deformation of the rock matrix, thus constructing the stress field equations. The fluid flow within the fracture is described based on the simplified Navier-Stokes equations, the flow field equations are constructed, and the relationship between fluid pressure and fracture volume change is established. Based on Fourier's law and the energy conservation equation, the heat exchange between fluids and solids is described, and the temperature field equation is constructed. Based on the stress field equation, flow field equation and temperature field equation, the thermal strain increment and pore pressure change caused by thermal expansion are introduced to construct a thermal-fluid-solid coupling mechanism, and the thermal stress increment and pore pressure increment are coupled into the mechanical constitutive equation to establish the thermal-fluid-solid-damage multiphysics field coupling control equation. Using Coulomb's law of friction as the basic model of fracture, a crack activation criterion is established. When the normal stress exceeds the tensile strength, it is judged as tensile failure; when the shear stress exceeds the maximum allowable shear force based on cohesion, internal friction angle and normal stress, it is judged as shear failure. The incremental normal and tangential stresses caused by thermal stress are superimposed on Coulomb's friction law to correct the critical condition for crack activation; once crack activation is determined, the tensile strength and cohesion are set to zero to simulate the damage evolution process. Based on mechanical and seepage properties, and using the multi-physics field coupling control equations of heat-fluid-solid-damage, a fracture activation model for fluid injection in fractured shale reservoirs is constructed. Numerical simulation of the fracture activation and slip process is performed to obtain multi-scale fracture activation and propagation slip maps during fracturing injection.

6. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 5, characterized in that, The thermo-fluid-solid-damage multiphysics coupling control equation is as follows: ； in, and It is the first step in the fluid injection process. i The crack to the first j The stress tensor components and strain tensor components corresponding to each crack; This is the Kronecker operator notation, and the table below represents the ij-th component of a symmetric tensor matrix, where... ; b It is the Biot coefficient; The shear modulus of the rock; The bulk modulus of the rock. , E The elastic modulus of the rock. v The Poisson's ratio of the rock; It is the linear thermal expansion coefficient of the rock; This is the sum of the strain tensors of the rock in different directions.

7. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 5, characterized in that, The simulation and characterization of the fracture activation and slip process also includes determining the fracture activation distribution characteristics based on the obtained multi-scale fracture activation and propagation slip maps during hydraulic fracturing injection, combined with microseismic event evaluation. Specifically, this includes: According to the clustering mechanism, when two slip segments in the multi-scale crack activation and extension slip map are adjacent to each other and the slip duration overlaps, the slip nodes within the preset radius R are merged and the microseismic slip segments are updated. The slip velocity is determined based on the updated microseismic slip segment. When the slip velocity is greater than the velocity threshold, it is determined to be a seismic slip; otherwise, it is determined to be a quiescent slip. For the slip segment determined to be seismically slipped, the seismic moment and magnitude are calculated based on the slip amount of the crack, the slip fracture area and the rock shear modulus, and microseismic spheres are generated. Based on the spatial location of microseismic spheres, a one-to-one correspondence between fault slip rate and microseismic events is established to determine the spatial distribution characteristics of fault activation.

8. The method for simulating the fluid injection fracture activation process in fractured shale reservoirs according to claim 1, characterized in that, The acquisition of well logging imaging data, microseismic data, and seismic data from shale reservoirs in the study area, and the construction of a multi-scale in-situ fracture parameter feature dataset, specifically includes: Collect well logging imaging data, microseismic data, and seismic data from oil and gas field construction sites; Grayscale and binary processing were performed on the logging imaging data to obtain in-situ characteristic parameters of well perimeter microfractures at the micrometer to centimeter level; Based on microseismic data and the in-situ characteristic parameters of meter-scale cracks were obtained using the inversion solution derived from the focal mechanism. The in-situ characteristic parameters of 100-meter-level fractures were obtained by processing seismic data using the ant body method. By fusing in-situ feature parameters of micro-fractures, micro-fractures, and fractures at the micrometer to centimeter scale, as well as those at the meter and hundred-meter scale, a multi-scale fracture in-situ parameter feature dataset is constructed.

9. A simulation system for the fluid injection fracture activation process in fractured shale reservoirs, characterized in that, include: A multi-scale dataset construction module is used to acquire well logging imaging data, microseismic data and seismic data of shale reservoirs in the study area, and to construct a multi-scale fracture in-situ parameter feature dataset. The parameter inversion module is used to input the multi-scale fracture in-situ parameter feature dataset into the pre-trained inversion model, invert the multi-scale fracture in-situ feature parameters, and obtain the inversion results; wherein, the inversion model is based on a physical information neural network, and a residual loss term coupled with shale seepage law and mechanical strength constraint is added after the output layer; The shale reservoir model construction module is used to construct a three-dimensional multi-scale fracture geometric model; integrate the inversion results into the geometric model to construct a three-dimensional multi-scale discrete fracture network model; construct a macroscopic heterogeneous shale reservoir matrix model based on the distribution function of macroscopic rock mechanical parameters of shale; and embed the three-dimensional multi-scale discrete fracture network model into the macroscopic heterogeneous shale reservoir matrix model to form a shale reservoir model that reflects the multi-scale discrete fracture network system and reservoir heterogeneity. The fracture activation and slip process simulation module is used to establish a fracture activation criterion based on the shale reservoir model and the Coulomb friction law. It also comprehensively considers the thermal stress, thermal expansion, in-situ shale reservoir environmental stress field, and the multi-physics field coupling effect of heat-fluid-solid-damage caused by external fluid injection, and constructs a fracture activation model of fractured shale reservoirs by fluid injection to simulate and characterize the fracture activation and slip process.

Images