A coupling simulation method for interwell interference and fracture network evolution in multi-well group time sequence fracturing process

Abstract

The present application belongs to the technical field of oil and gas field development, and particularly relates to a coupling simulation method for inter-well interference and fracture network evolution in multi-well group time sequence fracturing process. The method firstly constructs a high-precision three-dimensional geomechanical model through multi-source data fusion, and solves the initial stress field of the reservoir based on the finite element method; then, time sequence fracturing multi-physical field coupling simulation is carried out according to the actual construction sequence, and the stress field is calculated in real time based on the finite element method; the updated stress field is taken as the initial condition for subsequent well simulation, so as to realize accurate characterization of stress interference across time sequence; finally, through the whole well group iteration simulation process, the numerical simulation of all fracturing wells is completed, and the fracture network characteristic parameters and stress interference cumulative effect are output. The method innovatively realizes the whole-process coupling analysis from geologic modeling to whole-well group fracturing simulation, and provides a reliable technical means for optimizing fracturing design.

Description

Technical Field

[0001] This invention belongs to the technical field of oil and gas field development, specifically relating to a coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing. Background Technology

[0002] With the large-scale development of unconventional oil and gas resources, horizontal well multi-stage fracturing technology has become a core method for developing low-permeability reservoirs such as tight oil and shale gas. However, in actual development, multi-well fracturing faces key technical challenges such as severe inter-well interference and uncontrollable fracture network expansion. Traditional fracturing design methods are mainly based on static geological models and single-well simulations, which are difficult to accurately characterize the dynamic interference effects during multi-well time-series fracturing, leading to problems such as low reservoir stimulation efficiency and rapid production decline.

[0003] Currently, numerical simulation studies on multi-well fracturing mainly focus on two aspects: one is fracture propagation simulation based on geomechanical parameters, and the other is interference effect analysis based on production data. Existing methods generally have the following limitations: (1) They adopt the assumption of a fixed geostress field and fail to consider the dynamic evolution of geostress during fracturing construction and production; (2) They lack quantitative characterization of inter-well stress interference during time-series fracturing; (3) Fracture network propagation simulation and production dynamic analysis are disconnected, making it difficult to achieve full life cycle optimization design.

[0004] In reservoir characterization, conventional methods mainly rely on static geological parameters (such as rock mechanical properties and physical parameters) for classification and evaluation, failing to organically combine fracturing stimulation characteristics with dynamic production response. Although dynamic evaluation methods based on microseismic monitoring and production data analysis have been developed in recent years, the lack of effective multi-field coupling simulation means still makes it difficult to achieve synergistic optimization of fracturing parameters and development strategies.

[0005] In the prior art, patent CN120234922A proposed a method for quantifying the inter-well interference intensity of shale oil horizontal wells based on physical models. Although this method establishes an evaluation chart of interference intensity through large-scale three-wellbore physical model experiments and permeability tests, it still has the following significant shortcomings: First, its experimental design can only simulate simple fracturing scenarios with a limited number of wells, and cannot truly reflect the complex working conditions of time-series fracturing of multiple well groups in the field; second, this method only evaluates the interference effect through static permeability changes, and fails to fully consider the complex action mechanism of multi-field coupling such as dynamic evolution of geostress and multiphase flow in actual reservoirs.

[0006] Although patent CN120276066A innovatively proposes a three-dimensional well network interference evaluation method based on dynamic oil drainage volume, its research method still has obvious limitations: on the one hand, the method fails to establish an effective experimental-numerical simulation correlation method, making it difficult to accurately extrapolate physical experimental results to field applications; on the other hand, its evaluation model mainly relies on the single indicator of oil drainage volume and lacks accurate characterization of the three-dimensional expansion morphology of fracture network.

[0007] Patent CN119914231A discloses a physical simulation experimental method and apparatus for dealing with inter-well interference. However, the scale and complexity of the apparatus limit its ability to accurately reflect key physical processes such as competitive propagation of multiple fracture clusters and multiphase flow during production. More importantly, the evaluation system established by this method is not closely integrated with measurable engineering parameters in the field (such as far-field stress changes, microseismic events, and production dynamic data), making it difficult to directly guide fracturing design and interference control strategies in the field as an experimental method.

[0008] Patent CN118520996A discloses a method for quantitative characterization and risk prediction of inter-well interference in shale oil fracturing horizontal wells. Although it innovatively proposes an interference evaluation method based on RTA (Rate Transient Analysis) production decline analysis, its research method still has significant limitations: on the one hand, the method fails to establish an effective experimental-numerical simulation correlation method, making it difficult to accurately extrapolate physical experimental results to field applications; on the other hand, its evaluation model mainly relies on the single indicator of drainage volume, lacking a precise characterization of the three-dimensional propagation morphology of the fracture network. These methods all fail to solve the problem of dynamic simulation of geostress throughout the entire fracturing-to-production cycle, particularly lacking a quantitative description of key mechanisms such as the accumulation of stress shadow effects and competitive propagation of fracture networks during time-series fracturing.

[0009] Therefore, there is an urgent need to develop a new multi-well-group time-series fracturing simulation method that can achieve the following technical objectives: 1) accurate characterization of the dynamic evolution of the geostress field; 2) quantitative evaluation of the inter-well interference effect; 3) coupled analysis of fracture network expansion and production dynamics. Summary of the Invention

[0010] The purpose of this invention is to provide a coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, addressing the aforementioned deficiencies. By integrating multi-source data fusion, dynamic geostress modeling, and multi-physics field coupled simulation, this invention achieves a full-process coupled analysis from geological modeling to full-well-group fracturing simulation, providing a new technical means for the efficient development of unconventional oil and gas reservoirs.

[0011] The technical solution of this invention is as follows:

[0012] A coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing includes the following steps:

[0013] S1. Establish a geomechanical model by combining multi-source data.

[0014] S1.1 Based on field logging data, the logging curves are processed, and reservoir physical property data along the well trajectory are calculated. The reservoir physical property data includes clay content, porosity, fluid saturation, permeability, rock density, Young's modulus, and Poisson's ratio. The field logging data includes P-wave transit time (DTP, AC), S-wave transit time (DTS), density logging (DEN), neutron logging (CNL), natural gamma logging (GR), and resistivity logging (RT, RXO, ML).

[0015] S1.2. Import the three parameters of the well corresponding to the logging data in S1.1, namely the measurement depth (MD), inclination angle (INCL), and azimuth angle (AZIM), into the Petrel software to determine the well trajectory of the well corresponding to the logging data.

[0016] S1.3. Import the reservoir property data calculated in S1.1 into the Petrel software and place the reservoir property data on the corresponding well trajectory. Each reservoir property data has its own drilling number, and the data can be placed on the corresponding well trajectory by combining the corresponding well measurement depth.

[0017] Simultaneously, welltop data from different wells are imported onto the corresponding well trajectories. Welltop data is the measured depth of a well when it encounters different layers; in other words, it is the collection of the top depth and corresponding attribute information of each layer penetrated by the drill.

[0018] S1.4, Establish the construction surface:

[0019] First, by combining the well top data imported in S1.3, the positions of several points on the same plane are obtained; then, by interpolating from points to form a surface, the structural surface can be established.

[0020] As mentioned earlier, well point stratification data represents the measured depth of a well when it encounters different strata. Therefore, by combining the measured depths of multiple wells at the same stratum when they reach it during drilling, the positions of multiple points on that stratum can be obtained; finally, a structural surface can be established by interpolation using a point-to-surface approach.

[0021] The following three interpolation methods are commonly used: Minimum Curvature, suitable for sparse wells; Kriging, which considers spatial correlation and requires a variation function model; and Natural Neighbor, which is locally accurate and suitable for complex structures. The appropriate interpolation method can be selected based on the specific circumstances.

[0022] S1.5, Vertical layering:

[0023] Based on the structural surfaces obtained in S1.4, the number of grid divisions in the vertical direction is set; if the total number of structural surfaces is N, then (N-1) vertical layers are set, that is, the number of grid divisions in the vertical direction is N-1; where, the thickness of a single layer = the total thickness of N structural surfaces / the number of grid divisions in the vertical direction.

[0024] S1.6. Full reservoir interpolation is performed on the reservoir property data calculated in S1.1 based on Gaussian random function simulation (GRFS) to characterize the spatial variability of reservoir properties. Spatial correlation is described using a variogram, and a randomized realization conforming to geological laws is generated.

[0025] First, the reservoir property data calculated in S1.1 is discretized into a grid form using discretized well log upscale data, where the number of vertical layers of the grid is equal to the number of vertical grid divisions set in S1.5.

[0026] Then, based on the Petrophysical module, the reservoir property data calculated in S1.1 are interpolated to the full reservoir data to obtain the geomechanical model.

[0027] S2. Calculate and simulate the initial geostress field

[0028] First, calculate the vertical stress, minimum horizontal principal stress, and maximum horizontal principal stress respectively; as follows:

[0029] ① Vertical stress: Vertical stress mainly considers the pressure formed by the weight of the overlying rock strata.

[0030] ;

[0031] In the formula, σ ν —Vertical stress, MPa;

[0032] ρ —Rock density, kg / m³ 3 ;

[0033] g —Acceleration due to gravity, 9.81 m / s²2 .

[0034] ②Minimum horizontal principal stress: Calculated using Eaton's formula.

[0035] ;

[0036] In the formula, σ hmin —Minimum horizontal principal stress, MPa;

[0037] α —Biot coefficient, ranging from 0.6 to 1.0;

[0038] P p — Pore ​​pressure, MPa;

[0039] σ ν —Vertical stress, MPa;

[0040] ν —Poisson's ratio.

[0041] ③ Maximum horizontal principal stress: ;

[0042] In the formula, σ Hmax —Maximum horizontal principal stress, MPa;

[0043] σ hmin —Minimum horizontal principal stress, MPa;

[0044] K —Tectonic stress coefficient; corrected by fitting fracturing data;

[0045] σ ν —Vertical stress, MPa.

[0046] The structural stress coefficient K This is a key parameter characterizing the influence of regional tectonic stress on the horizontal geostress field. To improve the accuracy and geological adaptability of its value, a fitting and correction process based on in-situ fracturing test data is adopted, as follows:

[0047] First, based on regional geomechanical characteristics or engineering experience, the tectonic stress coefficient is determined. K Set reasonable initial values K 0, and based on this, calculate the theoretical maximum horizontal principal stress.

[0048] Subsequently, by combining hydraulic fracturing construction data from the target block or adjacent areas, the measured value of the maximum horizontal principal stress is obtained through pressure inversion or fracture monitoring.

[0049] Finally, based on the comparison of theoretical and measured values ​​across multiple wells and multiple sections, the tectonic stress coefficient was systematically adjusted.K This process allows the theoretically calculated values ​​to continuously approach the measured results, ultimately determining the optimal structural stress coefficient applicable to this block.

[0050] Then, the direction of the major axis of the wellbore collapse is observed through imaging logging (such as FMI), which is σ. hmin The direction of σ is perpendicular to the long axis of the well wall collapse. Hmax The direction.

[0051] Wellbore collapse is caused by the redistribution of stress around the wellbore after drilling, under the influence of the maximum horizontal principal stress (σ). Hmax In the direction of (i.e. parallel to the minimum horizontal principal stress σ) hmin The tangential stress concentration exceeds the shear strength of the rock, leading to symmetrical shear failure of the wellbore and the formation of an elliptical wellbore.

[0052] S3, Time-series fracturing multiphysics coupling simulation

[0053] S3.1 Fracturing fluid parameter settings:

[0054] Based on the type of fracturing fluid used in the on-site fracturing operation, the fracturing fluid parameters in the simulated work area are set. The fracturing fluid parameters include: fracturing fluid viscosity, frictional change of fracturing fluid with injection speed, and initial spurt loss of fracturing fluid.

[0055] Initial filtration loss characterizes the rapid filtration behavior of fracturing fluid when it comes into contact with the fracture wall and before the filter cake is fully formed, as well as the filtration loss coefficient of the fracturing fluid.

[0056] S3.2, Propionage parameter settings:

[0057] The proppant parameters include the proppant mesh size (particle size) and the type of proppant (e.g., ceramsite, silica sand, etc.). These parameters must correspond to the type of proppant used on site.

[0058] S3.3, Fracturing Pump Injection Program Settings:

[0059] The parameters of the fracturing pumping procedure include: the flow rate (pumping speed) for each construction stage, the type of fracturing fluid used, the total pumping volume, the type of proppant, and the proppant concentration.

[0060] The pumping procedure must strictly comply with the pumping procedure used during on-site construction.

[0061] S3.4 Simulation of crack propagation:

[0062] Based on the geomechanical model established in S1 and the fracturing program set in S3.1-3.3, a simulation of fracture network propagation was conducted using the Kinetix module in Petrel software. The fracture network propagation direction is along the direction of the maximum horizontal principal stress.

[0063] The pressure changes inside the fracture during the fracturing process were recorded using Petrel software, which was then used for subsequent simulation updates of the geostress field.

[0064] Petrel software calculates and analyzes the filtration loss of fracturing fluid based on its filtration properties, which is then used to update the saturation field.

[0065] Petrel software calculates the conductivity at different locations within the fracture based on the concentration and distribution of proppant within the fracture.

[0066] S3.5 Production simulation after fracturing:

[0067] If a horizontal well enters the flowback and production stage immediately after fracturing, meaning there is a construction process such as production of the mother well between the fracturing construction of the mother well and the fracturing construction of the daughter well, then the dynamic pressure field and seepage field formed during its production process will be further simulated.

[0068] First, the production grid is exported. The fracture morphology of the parent well based on the fracturing simulation results is exported to an unstructured grid. Petrel software calculates the reservoir physical parameters of the unstructured grid, including porosity and permeability.

[0069] Then, the Intersect simulator was used to simulate the flowback production of the parent well to obtain the reservoir pressure field at different production stages.

[0070] S3.6 Calculate and simulate the geostress field after fracturing / production.

[0071] Based on the finite element method, combined with the formation pressure field at different times, the changes in local geostress field caused by fracturing and production of the mother well are calculated in real time, and the stress field at different time points is calculated by the pressure field at different time points.

[0072] (1) Calculate the change of geostress field with respect to pressure field:

[0073] The geostress varying with pressure is calculated using the following formula:

[0074] ;

[0075] In the formula, σ yy —Stress component in the y direction;

[0076] σ xx —Stress component in the x direction;

[0077] σh0 —Initial minimum horizontal principal stress;

[0078] σ H0 —Initial maximum horizontal principal stress;

[0079] α —Biot coefficient, ranging from 0.6 to 1.0;

[0080] ν —Poisson's ratio;

[0081] P p — Changes in pore pressure.

[0082] The above formula is derived based on Biot's porosity elasticity equation and the strain-stress coupling relationship.

[0083] Pore ​​elastic constitutive relation:

[0084] ;

[0085] In the formula, — Tensor components of the current total stress;

[0086] — Tensor components in the initial state;

[0087] — Strain tensor components;

[0088] — Shear modulus describes a material's ability to resist shear deformation. ;

[0089] — Volumetric strain (satisfying: );

[0090] — Elastic property description parameters.

[0091] Expanding the above formula into component form:

[0092] ;

[0093] In the formula,

[0094] σ yy — y Stress components in the direction;

[0095] σxx — x Stress components in the direction;

[0096] σ zz — z Stress components in the direction;

[0097] G —Shear modulus, which describes a material's ability to resist shear deformation;

[0098] — x Strain components in the direction;

[0099] — y Strain components in the direction;

[0100] σ h0 —Initial minimum horizontal principal stress;

[0101] σ H0 —Initial maximum horizontal principal stress;

[0102] α —Biot coefficient, ranging from 0.6 to 1.0;

[0103] — Elastic property description parameters;

[0104] P p — Changes in pore pressure.

[0105] The vertical stress is mainly generated by the overlying rock strata, so its magnitude remains basically constant. Therefore, we can obtain the following relationship between volumetric strain and pore pressure:

[0106] ;

[0107] ;

[0108] In the formula, σ xx — x Stress components in the direction;

[0109] σ H0 —Initial maximum horizontal principal stress;

[0110] G —Shear modulus, which describes a material's ability to resist shear deformation;

[0111] — Elastic property description parameters;

[0112] — x Strain components in the direction;

[0113] — y Strain components in the direction;

[0114] α —Biot coefficient, ranging from 0.6 to 1.0;

[0115] P p — Pore ​​pressure variation;

[0116] σ yy — y Stress components in the direction;

[0117] σ h0 —Initial minimum horizontal principal stress.

[0118] Solving the two formulas simultaneously yields:

[0119] ;

[0120] In the formula, — x Strain components in the direction;

[0121] — y Strain components in the direction;

[0122] σ xx — x Stress components in the direction;

[0123] σ yy — y Stress components in the direction;

[0124] G —Shear modulus, which describes a material's ability to resist shear deformation;

[0125] σ H0 —Initial maximum horizontal principal stress;

[0126] σ h0 —Initial minimum horizontal principal stress.

[0127] Substituting both into the volumetric strain equation:

[0128] ;

[0129] In the formula, σ xx — xStress components in the direction;

[0130] σ H0 —Initial maximum horizontal principal stress;

[0131] G —Shear modulus, which describes a material's ability to resist shear deformation;

[0132] σ yy — y Stress components in the direction;

[0133] σ h0 —Initial minimum horizontal principal stress;

[0134] α —Biot coefficient, ranging from 0.6 to 1.0;

[0135] P p — Pore ​​pressure variation;

[0136] — Elastic property description parameters.

[0137] constant The expression:

[0138] ;

[0139] In the formula, — Elastic property description parameters;

[0140] G —Shear modulus, which describes a material's ability to resist shear deformation;

[0141] ν —Poisson's ratio.

[0142] Ultimately, we can obtain:

[0143] ;

[0144] σ xx — x Stress components in the direction;

[0145] σ yy — y Stress components in the direction;

[0146] σ H0 —Initial maximum horizontal principal stress;

[0147] σ h0 —Initial minimum horizontal principal stress;

[0148] α —Biot coefficient, ranging from 0.6 to 1.0;

[0149] ν —Poisson's ratio;

[0150] P p — Changes in pore pressure.

[0151] Assuming the boundary is a fixed displacement boundary, therefore ( ).

[0152] ;

[0153] In the formula, σ xx — x Stress components in the direction;

[0154] σ yy — y Stress components in the direction;

[0155] σ H0 —Initial maximum horizontal principal stress;

[0156] σ h0 —Initial minimum horizontal principal stress.

[0157] The final result is the change in the geostress field after the change in the pressure field:

[0158] ;

[0159] In the formula, σ yy — y Stress components in the direction;

[0160] σ xx — x Stress components in the direction;

[0161] σ h0 —Initial minimum horizontal principal stress;

[0162] σ H0 —Initial maximum horizontal principal stress;

[0163] α —Biot coefficient, ranging from 0.6 to 1.0;

[0164] ν —Poisson's ratio;

[0165] P p— Changes in pore pressure.

[0166] (2) Calculate the direction of geostress deflection

[0167] Calculate the principal stress direction angle using the following formula:

[0168] ;

[0169] ;

[0170] ;

[0171] In the formula, θ— Angle between the direction of maximum principal stress and the reference direction;

[0172] σ xx — x Stress components in the direction;

[0173] σ yy — y Stress components in the direction;

[0174] σ xy —Action perpendicular to x On the plane of the axis and parallel to the direction y The shear stress components of the shaft;

[0175] k — Stress variation non-uniformity coefficient;

[0176] Δσ xx — x Change in directional normal stress;

[0177] Δσ yy — y Change in directional normal stress;

[0178] P p — Changes in pore pressure.

[0179] By combining the calculation of the magnitude of geostress in (1) and the calculation of geostress deflection in (2), the connection between the pressure field and the geostress field is established, and the state of the geostress field at any time is finally obtained.

[0180] The above formula is derived as follows:

[0181] On a two-dimensional plane, the initial stress components are: If the coordinate system rotates by an angle Then the new stress components should be:

[0182] ;

[0183] In the formula, σ xx — x Stress components in the direction;

[0184] σ yy — y Stress components in the direction;

[0185] σ xy —Action perpendicular to x On the plane of the axis and parallel to the direction y The shear stress components of the shaft;

[0186] θ— The angle between the direction of the maximum principal stress and the reference direction.

[0187] The principal stress direction is the direction where the shear stress is zero, therefore .

[0188] We can obtain:

[0189] ;

[0190] ;

[0191] .

[0192] Combining the above formulas, we can obtain:

[0193] ;

[0194] In the formula, σ xx — x Stress components in the direction;

[0195] σ yy — y Stress components in the direction;

[0196] σ xy —Action perpendicular to x On the plane of the axis and parallel to the direction y The shear stress components of the shaft;

[0197] θ— The angle between the direction of the maximum principal stress and the reference direction.

[0198] Therefore, the new principal stress direction angle can be obtained as:

[0199] .

[0200] S4. Dynamic transfer of geostress field

[0201] First, the in-situ stress field calculated and simulated in S4 after fracturing / production is used as the initial condition for subsequent well simulations; the geomechanical data obtained from simulating the in-situ stress field after fracturing / production, including the maximum horizontal principal stress, minimum horizontal principal stress, vertical stress, and principal stress direction angle, are used as the fracturing reservoir conditions.

[0202] Then, using the "Promote to Kinetix zone set inputs" system, the initial geostress state is replaced with the latest geostress state. This ensures that the geostress state during subsequent well fracturing is the most recently simulated geostress state.

[0203] Thus, an algorithm for updating and transferring the geostress field was established, and a quantitative index system for stress disturbance was developed, including parameters such as the range of stress disturbance and the turning angle, to achieve accurate cross-time-series transfer of geostress evolution.

[0204] S5, Full-well-group iterative simulation

[0205] Following the construction sequence of the horizontal well group on site, repeat S3 and S4 cycles to conduct iterative simulation of the entire well group based on the process of "mother well fracturing / production - geostress field simulation - updating geostress field - daughter well fracturing" until all wells have been simulated, ensuring that the geostress field of each fracturing well is the latest geostress field at the time of fracturing.

[0206] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, wherein the clay content of reservoir physical property data in step S1.1 is calculated according to the following formula:

[0207] ;

[0208] In the formula, —Mud content;

[0209] —Gamma value of pure sandstone;

[0210] —Gamma value of pure mudstone;

[0211] —Well logging measurements.

[0212] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, wherein the porosity of reservoir physical property data in step S1.1 is calculated according to the following formula:

[0213] ;

[0214] In the formula, —Porosity;

[0215] —Measured acoustic time difference; ;

[0216] —Time difference of rock skeleton; For typical sandstone, the value is 55.5. The value for limestone is 47.5. The value for dolomite is 43.5. ;

[0217] —Mud content;

[0218] —Mudstone travel time;

[0219] —Fluid time difference; The value is usually 189. ;

[0220] — Compaction correction factor; where the value for undeveloped strata is 1.

[0221] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, in step S1.1, calculates the fluid saturation of reservoir physical property data according to the following formula (Archie formula):

[0222] ;

[0223] ;

[0224] In the formula, —Water saturation;

[0225] —hydrocarbon saturation;

[0226] —Lithology coefficient; value ranges from 0.6 to 1.0;

[0227] Φ m —Influencing factors related to porosity;

[0228] —Cementation index; value ranges from 1.8 to 2.2;

[0229] —Saturation index; typically taken as 2;

[0230] — Formation water resistivity ;

[0231] — Formation resistivity .

[0232] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, in step S1.1, calculates the permeability of reservoir physical property data according to the following formula:

[0233] ;

[0234] In the formula, —Penetration rate, mD;

[0235] —Porosity;

[0236] —Bound water saturation; can be obtained through core experiments.

[0237] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, wherein the rock density of reservoir physical property data in step S1.1 is calculated according to the following formula:

[0238] ;

[0239] In the formula, —Rock bulk density, g / cm³ 3 ;

[0240] —Electron density, e / cm 3 .

[0241] Density logging tools (such as LDT and FDC) emit medium-energy gamma rays (137Cs or 60Co source) into the formation. The rays undergo Compton scattering with electrons in the formation, and the detector records the intensity of the scattered gamma rays. There is a linear relationship between electron density and volume density.

[0242] The Young's modulus of the reservoir properties in step S1.1 is calculated according to the following formula:

[0243] ;

[0244] In the formula, —Dynamic Young's modulus, GPa;

[0245] —Rock bulk density, g / cm³ 3 ;

[0246] —P-wave velocity, m / s;

[0247] — Shear wave velocity, m / s.

[0248] The Poisson's ratio of the reservoir property data in step S1.1 is calculated according to the following formula:

[0249] ;

[0250] —Dynamic Poisson's ratio;

[0251] —P-wave velocity, m / s;

[0252] — Shear wave velocity, m / s.

[0253] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, in step S2, if there are already fracturing wells in the reservoir, the maximum and minimum horizontal principal stresses are obtained through the construction pressure:

[0254] During construction, the rupture pressure satisfies the Hubbert-Willis equation:

[0255] ;

[0256] In the formula, —The burst pressure (peak pressure) during construction can be obtained from the construction pressure curve, in MPa;

[0257] σ hmin —Minimum horizontal principal stress, MPa;

[0258] σ Hmax —Maximum horizontal principal stress, MPa;

[0259] P p — Pore ​​pressure, MPa;

[0260] —Tensile strength of rock, MPa; can be obtained from core experiments.

[0261] The mechanical equilibrium condition for crack closure is: ;

[0262] —The crack closure pressure during construction can be obtained from the construction pressure curve analysis;

[0263] σ hmin —Minimum horizontal principal stress, MPa.

[0264] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, wherein the filtrate loss in step S3 is calculated according to the following formula:

[0265] ;

[0266] In the formula, —Filtration loss per unit area ;

[0267] —Initial filtration loss per unit area ;

[0268] — Filtration coefficient ;

[0269] —Construction time, in minutes.

[0270] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, wherein the proppant concentration in step S3 is calculated using the following formula:

[0271] ;

[0272] —Propionate concentration, kg / m³ 3 ;

[0273] ρ p —Propion density, kg / m³ 3 .

[0274] During on-site construction, the proppant concentration is generally not given; in most cases, only the sand ratio is given. Therefore, the sand ratio and proppant concentration need to be converted according to the above formula.

[0275] In this invention, the coupled simulation method for inter-well interference and fracture network evolution during multi-well time-series fracturing can employ one of the following three commonly used interpolation methods in step S1.4: Minimum Curvature interpolation, suitable for sparse wells; Kriging interpolation, which considers spatial correlation and requires a variogram model; and Natural Neighbor interpolation, which is locally accurate and suitable for complex structures. The appropriate interpolation method can be selected based on the specific circumstances.

[0276] The beneficial effects of this invention are as follows: This invention provides a coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing. This method first constructs a high-precision three-dimensional geomechanical model through multi-source data fusion and solves the initial reservoir stress field based on the finite element method. Then, it conducts multi-physics coupled simulation of time-series fracturing according to the actual construction sequence, including fluid-structure interaction simulation of the fracturing process (considering characteristics such as fracturing fluid loss and proppant migration), dynamic pressure and seepage field simulation during the production stage (if a production stage exists), and real-time updating calculation of the stress field based on the finite element method. Furthermore, it establishes a dynamic stress field transfer mechanism, using the updated stress field as the initial condition for subsequent well simulations, achieving accurate cross-time characterization of stress interference. Finally, through a full-well group iterative simulation process, it completes the numerical simulation of all fracturing wells, outputting fracture network characteristic parameters and the cumulative effect of stress interference. This method innovatively realizes a full-process coupled analysis from geological modeling to full-well group fracturing simulation, providing a reliable technical means for optimizing fracturing design. Attached Figure Description

[0277] Figure 1 This is a flowchart illustrating the coupled simulation method of inter-well interference and fracture network evolution during multi-well-group time-series fracturing as described in this invention.

[0278] Figure 2 This is the well trajectory of the exploration well and the imported logging curve in Embodiment 1 of the present invention.

[0279] Figure 3 This is a layered information diagram of the well trajectory and its corresponding layer in Embodiment 1 of the present invention.

[0280] Figure 4 This is a layer diagram obtained based on layered information in Embodiment 1 of the present invention.

[0281] Figure 5 This is a detailed view of the vertical layering of the model in Embodiment 1 of the present invention.

[0282] Figure 6 The well logging curves in Embodiment 1 of the present invention are discretized and exported to a grid diagram.

[0283] Figure 7 This is a reservoir stress distribution diagram in Embodiment 1 of the present invention.

[0284] Figure 8 This is a schematic diagram of the wellbore under horizontal stress control in Embodiment 1 of the present invention.

[0285] Figure 9 This is a graph showing the frictional properties of the fracturing fluid used in Example 1 of the present invention.

[0286] Figure 10 This is a graph showing the filtration properties of the fracturing fluid in Embodiment 1 of the present invention.

[0287] Figure 11 This is a permeability diagram of the 20 / 40 mesh proppant in Example 1 of the present invention under different closure stresses.

[0288] Figure 12 This is a diagram showing the pressure distribution of fracturing fluid during fracturing operations in Embodiment 1 of the present invention.

[0289] Figure 13 This is a fracture network diagram obtained from the fracturing simulation in Embodiment 1 of the present invention.

[0290] Figure 14 This is an example of deriving an unstructured mesh diagram based on the crack network morphology in Embodiment 1 of the present invention.

[0291] Figure 15 This is a reservoir pressure field diagram after flowback production simulation in Embodiment 1 of the present invention.

[0292] Figure 16 This is a diagram showing the stress increase near the wellbore after fracturing in Embodiment 1 of the present invention.

[0293] Figure 17 This diagram illustrates the increase in geostress and the deflection of geostress direction following fracturing in Embodiment 1 of the present invention.

[0294] Figure 18 This is a diagram showing the fracturing results of well Le61-H711 in Embodiment 1 of the present invention without considering stress field changes.

[0295] Figure 19 The diagram shows the fracturing results of well Le61-H711 in Embodiment 1 of the present invention when considering stress field changes.

[0296] Figure 20 This is a fracture morphology diagram of well Le13-H704 in Embodiment 1 of the present invention without considering changes in the stress field.

[0297] Figure 21 The image shows the fracture morphology of well Le13-H704 in Embodiment 1 of the present invention after considering the change in stress field.

[0298] Figure 22 The diagram shows the mesh morphology of Le13-H707 in Embodiment 1 of the present invention when considering changes in the geostress field.

[0299] Figure 23 This is a diagram of the mesh morphology of Le13-H707 in Embodiment 1 of the present invention without considering changes in the geostress field. Detailed Implementation

[0300] The technical solution of the present invention will now be described in detail with reference to the accompanying drawings.

[0301] Example 1

[0302] The simulation method described in this invention is used to perform coupled simulation of inter-well interference and fracture network evolution during multi-well-group time-series fracturing. The specific steps are as follows:

[0303] S1. Establish a geomechanical model by combining multi-source data:

[0304] S1.1 Based on field logging data, the logging curves are processed, and reservoir physical property data along the well trajectory are calculated. The reservoir physical property data includes clay content, porosity, fluid saturation, permeability, rock density, Young's modulus, and Poisson's ratio. The field logging data includes P-wave transit time (DTP, AC), S-wave transit time (DTS), density logging (DEN), neutron logging (CNL), natural gamma logging (GR), and resistivity logging (RT, RXO, ML).

[0305] Actual well logging involves a series of data. The calculation process will now be explained using a specific logging point as an example:

[0306] (1) Calculation of clay content:

[0307] .

[0308] (2) Porosity calculation:

[0309] .

[0310] (3) Calculation of fluid saturation:

[0311] ① Water saturation: .

[0312] ② Hydrocarbon saturation: .

[0313] (4) Permeability calculation: .

[0314] (5) Rock density calculation: .

[0315] (6) Calculation of Young's modulus:

[0316] .

[0317] (7) Calculation of Poisson's ratio: .

[0318] S1.2. Import the three parameters of the well corresponding to the logging data in S1.1, namely the measurement depth (MD), inclination angle (INCL), and azimuth angle (AZIM), into the Petrel software to determine the well trajectory of the well corresponding to the logging data.

[0319] Taking the well trajectory of a vertical well (pilot well) as an example, the data in Table 1 is imported into the Petrel software in the form of a txt text file.

[0320] Table 1

[0321]

[0322] Note: This diameter represents a perfectly vertical exploration well, therefore its inclination angle and inclination azimuth are both 0.

[0323] S1.3. Import the reservoir property data calculated in S1.1 into the Petrel software and place the reservoir property data on the corresponding well trajectory. Each reservoir property data has its own drilling number, and the data can be placed on the corresponding well trajectory by combining the corresponding well measurement depth.

[0324] Simultaneously, welltop data from different wells are imported onto the corresponding well trajectories. Welltop data is the measured depth of a well when it encounters different layers; in other words, it is the collection of the top depth and corresponding attribute information of each layer penetrated by the drill.

[0325] Based on the aforementioned logging data, the logging curves are processed, and the obtained reservoir physical property data (porosity, permeability, saturation, etc.) are imported into the Petrel software and placed on the corresponding well trajectory according to the corresponding well number.

[0326] Taking water saturation as an example, the well logging interpretation results are as follows: Figure 2 As shown.

[0327] Along the well trajectory, wellpoint stratification data from different wells are imported into the well trajectory. Currently, only two layers of stratification information have been imported, meaning that wellpoint stratification data exists on the well trajectory of each well. For example... Figure 3 As shown, there are two points on each well trajectory, which are the specific locations where the two layers pass through the well trajectory.

[0328] S1.4, Establish the construction surface:

[0329] First, by combining the well point layering data imported in S1.3, the positions of several points on the same layer are obtained; then, by interpolating from points to form a surface, the structural surface can be established.

[0330] As mentioned earlier, well point stratification data represents the measured depth of a well when it encounters different strata. Therefore, by combining the measured depths of multiple wells at the same stratum when they reach it during drilling, the positions of multiple points on that stratum can be obtained; finally, a structural surface can be established by interpolation using a point-to-surface approach.

[0331] Because the pilot wells are relatively sparse, the Minimum Curvature (a smooth surface suitable for sparse wells) interpolation method was chosen for this study. Figure 4 As shown.

[0332] S1.5, Vertical layering:

[0333] Based on the two structural surfaces obtained in S1.4, one vertical layer is required this time. The thickness of a single layer is: total thickness / number of grid divisions. Currently, the reservoir thickness is 20m, therefore 20 layers are set, resulting in an average thickness of 1m per layer. Figure 5 As shown.

[0334] S1.6. Perform full reservoir interpolation on the reservoir property data calculated in S1.1 based on Gaussian random function simulation (GRFS).

[0335] First, Well log upscaling is used to discretize the calculated reservoir property data into a grid. The vertical number of grid layers is equal to the 20 layers set earlier. For example... Figure 6 As shown.

[0336] Then, based on the Petrophysical module, the calculated reservoir property data is interpolated to the data for the entire block. Taking the distribution of the maximum horizontal principal stress in the reservoir as an example, such as... Figure 7 As shown.

[0337] S2. Calculate and simulate the initial geostress field:

[0338] First, calculate the vertical stress, minimum horizontal principal stress, and maximum horizontal principal stress respectively; as follows:

[0339] ① Vertical stress: The density of the overlying rock strata is simplified to an average of 2900 kg / m³. 3 .

[0340] .

[0341] ② Minimum horizontal principal stress: Given that the target reservoir pore pressure is 30.3 MPa, the formation fracturing pressure based on historical fracturing data from adjacent blocks is 60.35 MPa, and the rock tensile strength measured by core experiments is 3.5 MPa.

[0342] .

[0343] ③ Maximum horizontal principal stress:

[0344] .

[0345] The direction of the long axis of the collapse observed by imaging logging (FMI) is σ.hmin The direction of , and the direction perpendicular to it is σ. hmax The direction.

[0346] Therefore, based on the on-site wellbore diagram, such as Figure 8 As shown, the direction of the maximum principal stress is close to 75° east of north.

[0347] S3, Time-series fracturing multiphysics coupling simulation

[0348] S3.1 Fracturing fluid parameter settings:

[0349] Based on the on-site fracturing procedure report, slick water was primarily used as the fracturing fluid in this operation. The viscosity and filtration loss coefficient of the fracturing fluid are shown in Table 2.

[0350] Table 2

[0351]

[0352] The change in friction of fracturing fluid with injection rate, such as Figure 9 As shown.

[0353] Initial filtrate loss of fracturing fluid, such as Figure 10 As shown.

[0354] S3.2, Propionage parameter settings:

[0355] According to the on-site fracturing report, 20 / 40 mesh natural sand was mainly used as proppant in this fracturing operation. The permeability of 20 / 40 mesh proppant under different closure stresses is as follows: Figure 11 As shown. The permeability of the proppant under different closure stresses is taken into account when deriving the unstructured mesh. Based on the closure stress at different locations of the crack, the permeability corresponding to the unstructured mesh can be obtained according to the "closure stress-permeability" relationship curve of the proppant.

[0356] S3.3, Fracturing Pump Injection Program Settings:

[0357] In on-site construction, proppant concentration is generally not specified; in most cases, only the sand ratio is given. Therefore, it is necessary to convert between the sand ratio and the proppant concentration. The conversion formula is as follows:

[0358] .

[0359] The on-site fracturing procedure was converted into a format that Petrel can recognize, as shown in Table 3.

[0360] Table 3

[0361]

[0362] After setting up the fracturing program, you can start the fracturing simulation of the first fracturing well.

[0363] S3.4 Simulation of crack propagation:

[0364] Based on the previously established geomechanical model and the pre-configured fracturing program, a simulation of fracture network propagation was conducted using the Kinetix module in Petrel software. The fracture network propagated along the direction of the maximum horizontal principal stress.

[0365] The pressure changes inside the fracture during the fracturing process were recorded using Petrel software, which was then used for subsequent simulation updates of the geostress field.

[0366] Petrel software calculates and analyzes the filtration loss of fracturing fluid based on its filtration properties, which is then used to update the saturation field.

[0367] Petrel software calculates the conductivity at different locations within the fracture based on the concentration and distribution of proppant within the fracture.

[0368] The pressure distribution within the fracture network is as follows: Figure 12 As shown in Table 4.

[0369] Depend on Figure 10 It is known that the initial filtration loss per unit area of ​​the fracturing fluid under the target reservoir is 0.00035, and the simulation results show that the fracture wall area is 2650 m². 2 The calculation method for fluid loss during fracturing is as follows:

[0370]

[0371] The total fracturing fluid usage can be obtained by summing the fluid volume columns in Table 3: 1377.61 m³. 3 .

[0372] Fracturing fluid filtrate ratio:

[0373]

[0374] The Petrel software can automatically identify and calculate the average pressure inside the fracture and the fracture conductivity, with the fracture conductivity calculation referenced in the figure. Figure 11 The software automatically based on Figure 11 The permeability distribution inside the crack is calculated using the curve pattern.

[0375] The software simulation results show that the average width of the fracture is 2.45 mm. Simultaneously, based on the "different closure stress-permeability" relationship curves we provided, the software calculates the average permeability within the entire fracture to be 185320 mD. Therefore, the average conductivity of the fracture can be calculated.

[0376] Table 4

[0377]

[0378] S3.5 Production simulation after fracturing:

[0379] First, the production grid is exported by extracting the fractures from the parent well into an unstructured grid. Petrel software is used to calculate reservoir properties, including porosity and permeability, from this unstructured grid. The post-fracturing fracture network is shown below. Figure 13 Export to unstructured production grids such as Figure 14 As shown.

[0380] Then, the Intersect simulator was used to simulate the flowback production of the parent well, obtaining the reservoir pressure field at different production stages. For example... Figure 15 As shown.

[0381] S3.6 Calculate and simulate the geostress field after fracturing / production.

[0382] (1) Calculate the change of geostress field with respect to pressure field:

[0383] After fracturing, the average pressure within the fracture reached 50.13 MPa. The stress variation near the fracture can then be calculated using Biot's porosity elasticity equation and the strain-stress coupling relationship.

[0384] .

[0385] After production, the pressure at a certain point decreased by 10 MPa, while the original formation pressure was 30.3 MPa. The stress change at that point can then be calculated using Biot's porosity elasticity equation and the strain-stress coupling relationship.

[0386] .

[0387] Fracturing causes an increase in stress near the wellbore, such as Figure 16 As shown.

[0388] (2) Calculate the direction of geostress deflection

[0389] During fracturing, if the maximum horizontal principal stress at a certain point is 50 MPa and the minimum horizontal principal stress is 48 MPa, then the deflection angle is (the negative sign here only represents the direction of the force, and the negative sign means opposite to the positive direction of the coordinate axis):

[0390]

[0391]

[0392] .

[0393] After fracturing, the ground stress is increased and the direction of the ground stress is deflected (arrow direction), such as... Figure 17 As shown.

[0394] S4. Dynamic transfer of geostress field

[0395] First, the in-situ stress field calculated and simulated in S4 after fracturing / production is used as the initial condition for subsequent well simulations; the geomechanical data obtained from simulating the in-situ stress field after fracturing / production, including the maximum horizontal principal stress, minimum horizontal principal stress, vertical stress, and principal stress direction angle, are used as the fracturing reservoir conditions.

[0396] Then, using the “Promote to Kinetix zone set inputs” system, the initial geostress state is replaced with the latest geostress state.

[0397] S5, Full-well-group iterative simulation

[0398] Following the construction sequence of the horizontal well group on site, repeat S3 and S4 cycles to conduct iterative simulation of the entire well group based on the process of "mother well fracturing / production - geostress field simulation - updating geostress field - daughter well fracturing" until all wells have been simulated, ensuring that the geostress field of each fracturing well is the latest geostress field at the time of fracturing.

[0399] After considering the change in stress field, the reservoir stimulation volume (SRV) of multiple wells changed significantly, as shown in Table 5.

[0400] Table 5

[0401]

[0402] Among them, Le61-H710 and Le61-H712 were the first to be fracturing, and since they are far apart, both were simulated on the original geostress field and did not affect each other. Therefore, the reservoir stimulation volume before and after them did not change.

[0403] Well Le61-H711 is located between H710 and H712, and its construction time was later than both of them. Therefore, the high-pressure zone formed by fracturing the two parent wells first led to the formation of high-stress zones on both sides. The high-stress zone hindered the propagation of the fracture in H711, thus reducing its average fracture length from 165.12m to 140.20m. Wells H702 and H701 exhibited the same situation for the same reason.

[0404] Without considering stress field variations, the average seam length is 165.12m. Figure 18 As shown.

[0405] Considering the stress field variation, the average seam length is 140.20m. Figure 19 As shown.

[0406] Because well Le13-H704 overlaps with the older well H712 in the vertical direction, and the stress angle of the fractured well H712 was significantly deflected after fracturing, the fracture propagation direction of the overlapping section was also significantly deflected. Well H703 also exhibited the same deflection phenomenon for the same reason.

[0407] All seams are straight, without considering changes in the stress field, such as Figure 20 As shown.

[0408] Considering the significant deflection of the crack after the stress field changes, such as Figure 21 As shown.

[0409] There are producing fracturing wells near Le13-H707, and the production caused a decrease in pressure, which in turn created a low-stress zone. This resulted in a reduction in the complexity of the fracturing fractures in well H707, and a more uneven extension on both sides.

[0410] Considering the reduced complexity after the change in the geostress field, such as Figure 22 As shown.

[0411] Without considering changes in the geostress field, the complexity of the seam mesh is high, such as... Figure 23 As shown.

Claims

1. A coupled simulation method for inter-well interference and fracture network evolution during multi-well-group time-series fracturing, characterized in that, Includes the following steps: S1. Establish a geomechanical model by combining multi-source data. S1.1 Based on the field logging data, the logging curves are processed, and reservoir physical property data along the well trajectory are calculated; among which, the reservoir physical property data include clay content, porosity, fluid saturation, permeability, rock density, Young's modulus, and Poisson's ratio; S1.

2. Import the three parameters of the well corresponding to the logging data in S1.1, namely the measurement depth, well inclination angle and well inclination azimuth angle, into the Petrel software to determine the well trajectory of the well corresponding to the logging data. S1.

3. Import the reservoir property data calculated in S1.1 into the Petrel software and place the reservoir property data on the corresponding well trajectory; At the same time, the well point stratification data of different wells are imported into the corresponding well trajectories; S1.4, Establish the construction surface: First, by combining the well point layering data imported in S1.3, the positions of several points on the same layer are obtained; then, by interpolating from points to surfaces, the structural surface can be established. S1.5, Vertical layering: Based on the structural surfaces obtained in S1.4, the number of grid divisions in the vertical direction is set; if the total number of structural surfaces is N, then (N-1) vertical layers are set, that is, the number of grid divisions in the vertical direction is (N-1); where, the thickness of a single layer = the total thickness of N structural surfaces / the number of grid divisions in the vertical direction; S1.

6. Perform full reservoir interpolation on the reservoir property data calculated in S1.1 based on Gaussian random function simulation; First, the continuous reservoir property data calculated in S1.1 is discretized into a grid by using discretized logging curve data. The number of vertical layers of the grid is equal to the number of vertical grid divisions set in S1.

5. Then, based on the rock physics module, the reservoir physical property data calculated in S1.1 is interpolated into full reservoir data to calculate the geomechanical model; S2. Calculate and simulate the initial geostress field First, calculate the vertical stress, minimum horizontal principal stress, and maximum horizontal principal stress respectively; as follows: ① Vertical stress: ; In the formula, σ ν —Vertical stress, MPa; ρ— Rock density, kg / m³ 3 ; g —Acceleration due to gravity, 9.81 m / s² 2 ; ②Minimum horizontal principal stress: ; In the formula, σ hmin —Minimum horizontal principal stress, MPa; α —Biot coefficient, ranging from 0.6 to 1.0; P p — Pore ​​pressure, MPa; σ ν —Vertical stress, MPa; ν —Poisson's ratio; ③ Maximum horizontal principal stress: ; In the formula, σ Hmax —Maximum horizontal principal stress, MPa; σ hmin —Minimum horizontal principal stress, MPa; K —Structural stress coefficient; σ ν —Vertical stress, MPa; Then, the direction of the major axis of the wellbore collapse, denoted as σ, is observed through imaging logging. hmin The direction of σ is perpendicular to the long axis of the well wall collapse. Hmax The direction; S3, Time-series fracturing multiphysics coupling simulation S3.1 Fracturing fluid parameter settings: Based on the type of fracturing fluid used in the on-site fracturing operation, the fracturing fluid parameters in the simulated work area are set; the fracturing fluid parameters include: fracturing fluid viscosity, the change in friction of the fracturing fluid with the injection rate, and the initial filtration loss of the fracturing fluid. S3.2, Propionage parameter settings: The proppant parameters include the mesh size and type of proppant; S3.3, Fracturing Pump Injection Program Settings: The parameters of the fracturing pump injection procedure include: the flow rate at each construction stage, the type of fracturing fluid used, the total injection volume, the type of proppant, and the proppant concentration; S3.4 Simulation of crack propagation: Combining the geomechanical model established in S1 and the fracturing program set in S3.1-3.3, a simulation of the propagation of the fracture network was carried out based on the Kinetix module in Petrel software; the fracture network propagation direction is along the direction of the maximum horizontal principal stress. The pressure changes inside the fracture during the fracturing process were recorded using Petrel software, which was then used for subsequent simulation updates of the geostress field. Petrel software calculates and analyzes the filtration loss of fracturing fluid based on its filtration properties, which is then used to update the saturation field. Petrel software calculates the conductivity at different locations inside the fracture based on the concentration and distribution of proppant inside the fracture. S3.5 Production simulation after fracturing: If a horizontal well enters the flowback and production stage immediately after fracturing, the dynamic pressure field and seepage field formed during its production process will be further simulated. First, the production grid is exported. The fracture morphology of the parent well is exported to the unstructured grid based on the fracturing simulation results. The Petrel software calculates the reservoir physical parameters of the unstructured grid, including porosity and permeability. Then, the Intersect reservoir numerical simulator was used to simulate the flowback production of the parent well and obtain the reservoir pressure field at different production stages. S3.6 Calculate and simulate the in-situ stress field after fracturing / production: Based on the finite element method, combined with the formation pressure field at different times, the changes in local geostress field caused by fracturing and production of the parent well are calculated in real time, and the stress field at different time points is calculated by the pressure field at different time points. (1) Calculate the change of geostress field with respect to pressure field: The geostress varying with pressure is calculated using the following formula: In the formula, σ yy —Stress component in the y direction; σ xx —Stress component in the x direction; σ h0 —Initial minimum horizontal principal stress; σ H0 —Initial maximum horizontal principal stress; α —Biot coefficient, ranging from 0.6 to 1.0; ν —Poisson's ratio; P p — Pore ​​pressure change, MPa; (2) Calculate the direction of geostress deflection: Calculate the principal stress direction angle using the following formula: ； ； ； In the formula, θ— Angle between the direction of maximum principal stress and the reference direction; σ xx — x Stress components in the direction; σ yy — y Stress components in the direction; σ xy —Action perpendicular to x On the plane of the axis and parallel to the direction y The shear stress components of the shaft; k — Stress variation non-uniformity coefficient; Δσ xx — x Change in directional normal stress; Δσ yy — y Change in directional normal stress; P p — Pore ​​pressure change, MPa; By combining the calculation of the magnitude of geostress in (1) and the calculation of geostress deflection in (2), the connection between the pressure field and the geostress field is established, and the state of the geostress field at any time is finally obtained. S4. Dynamic transfer of geostress field First, the in-situ stress field calculated and simulated in S4 after fracturing / production is used as the initial condition for subsequent sub-well simulations; the geomechanical data obtained from simulating the in-situ stress field after fracturing / production, including the maximum horizontal principal stress, minimum horizontal principal stress, vertical stress, and principal stress direction angle, is used as the fracturing reservoir condition. Then, the initial geostress state is replaced with the latest geostress state using the Promote to Kinetix zone set inputs system. S5, Full-well-group iterative simulation Following the construction sequence of the horizontal well group on site, repeat S3 and S4 cycles to conduct iterative simulation of the entire well group based on the process of "mother well fracturing / production - geostress field simulation - updating geostress field - daughter well fracturing" until all wells have been simulated, ensuring that the geostress field of each fracturing well is the latest geostress field at the time of fracturing.

2. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing as described in claim 1, characterized in that, The clay content of the reservoir physical property data in step S1.1 is calculated according to the following formula: ； In the formula, —Mud content; —Gamma value of pure sandstone; —Gamma value of pure mudstone; —Well logging measurements.

3. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The porosity of the reservoir properties in step S1.1 is calculated according to the following formula: ； In the formula, —Porosity; —Measured acoustic time difference; ; —Rock skeleton time difference; ; —Mud content; —Mudstone travel time; —Fluid time difference; ; —Compaction correction coefficient.

4. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The fluid saturation of the reservoir property data in step S1.1 is calculated according to the following formula: ； ； In the formula, —Water saturation; —hydrocarbon saturation; —Lithology coefficient; Φ m —Influencing factors related to porosity; —Cementation index; —Saturation index; — Formation water resistivity ; — Formation resistivity .

5. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The permeability of the reservoir properties data in step S1.1 is calculated according to the following formula: ； In the formula, —Penetration rate, mD; —Porosity; —Bound water saturation.

6. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The rock density of the reservoir physical property data in step S1.1 is calculated according to the following formula: ； In the formula, —Rock bulk density, g / cm³ 3 ; —Electron density, e / cm 3 ; The Young's modulus of the reservoir properties in step S1.1 is calculated according to the following formula: In the formula, —Dynamic Young's modulus, GPa; —Rock bulk density, g / cm³ 3 ; —P-wave velocity, m / s; —Shear wave velocity, m / s; The Poisson's ratio of the reservoir property data in step S1.1 is calculated according to the following formula: ； In the formula, —Dynamic Poisson's ratio; —P-wave velocity, m / s; — Shear wave velocity, m / s.

7. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, In step S2, if there are already fractured wells in the reservoir, the maximum and minimum horizontal principal stresses are obtained through the construction pressure: During construction, the rupture pressure satisfies the Hubbert-Willis equation: ； In the formula, —Rupture pressure during construction, MPa; σ hmin —Minimum horizontal principal stress, MPa; σ Hmax —Maximum horizontal principal stress, MPa; P p — Pore ​​pressure, MPa; —Tensile strength of rock, MPa; The mechanical equilibrium condition for crack closure is: ; In the formula, —Crack closure pressure during construction; σ hmin —Minimum horizontal principal stress, MPa.

8. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The filtration loss per unit area of ​​the filtrate in step S3 is calculated according to the following formula: ； In the formula, —Filtration loss per unit area ; —Initial filtration loss per unit area ; — Filtration coefficient ; —Construction time, in minutes.

9. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, The proppant concentration in step S3 is calculated using the following formula: ； In the formula, —Propionate concentration, kg / m³ 3 ; ρ p —Propion density, kg / m³ 3 .

10. The coupled simulation method for inter-well interference and fracture network evolution during multi-well group time-series fracturing according to claim 1, characterized in that, In step S1.4, the interpolation method for establishing the construction surface is one of smooth surface, geospatial interpolation, or natural neighbor interpolation.

Images