CN113919196A - Reservoir three-dimensional stress field simulation method, simulation system, terminal and storage medium - Google Patents

Abstract

The invention belongs to the technical field of oil and gas field development geology and petroleum engineering rock mechanics, and discloses a reservoir three-dimensional stress field simulation method, a reservoir three-dimensional stress field simulation system, a reservoir three-dimensional stress field simulation terminal and a reservoir three-dimensional stress field storage medium. And classifying the seismic facies according to the development characteristics of the bedding fractures and the natural fractures on the well. Under the constraint of the seismic facies, the crack lines are tracked by utilizing the enhanced curvature property, a discrete natural crack network deterministic model is established on the seismic scale, meanwhile, the influences of local strength weakening and stress concentration caused by geological discontinuous surfaces such as faults and natural cracks on the whole rock mechanical field and the ground stress field are considered, a related grid conversion program is developed to realize the grid conversion of the geological model and the ground stress finite element model, and the optimum solution acquisition of the ground stress field is realized by utilizing a competitive particle swarm optimization algorithm. The method obtains the three-dimensional ground stress result with high precision and resolution, and provides important technical support for subsequent oil and gas field exploration and development work such as engineering dessert prediction, well structure design of drilling, fracturing yield-increasing scheme analysis and the like.

Description

Reservoir three-dimensional stress field simulation method, simulation system, terminal and storage medium

Technical Field

The invention belongs to the technical field of oil and gas field development geology and petroleum engineering rock mechanics, and particularly relates to a reservoir three-dimensional stress field simulation method, a reservoir three-dimensional stress field simulation system, a reservoir three-dimensional stress field simulation terminal and a reservoir three-dimensional stress field storage medium.

Background

The crustal stress is mainly formed by superposition of self-weight stress and structural stress, is a fundamental acting force causing deformation and destruction of underground engineering, and has important significance on series of problems such as reservoir dessert evaluation, well track design, well wall instability, sand production, hydraulic fracturing design, casing damage, well pattern optimization and the like by accurately mastering the size and the direction of the crustal stress in exploration and development of petroleum and natural gas.

The method has the following problems in the aspects of setting and simulating the three-dimensional ground stress field that firstly, the influence of geological structures such as faults, natural cracks and the like on the stress field is less considered; secondly, the description of the rock mechanics field and the natural fracture mostly depends on the inter-well interpolation calculation, and the accuracy of the inter-well parameters is deficient; the geological grid is not matched with the finite element grid, parameters cannot be effectively transmitted, and grid distortion is easy to occur in complex geological regions such as stratum interfaces, fault planes, pinch-out and the like; and fourthly, the ground stress field with correct distribution rule can be simulated, but the traditional algorithms such as iterative regression and the like adopted in the process are easy to fall into the local optimal solution, and the convergence speed and the inversion precision effect are general.

The geological or oil reservoir model is mostly modeled by adopting an angular point grid, has stronger flexibility and can better represent the heterogeneity of the oil-gas reservoir. The finite element method has great advantages in solving nonlinear mechanics problems and complex stress strain problems, and the stress distribution condition of the target area can be accurately obtained through finite element ground stress simulation based on the finite element grid. However, the geological, reservoir and finite element models are usually not matched with each other, and the mutual identification and conversion accuracy is poor, so that the simplification of the process is required when the ground stress finite element simulation is performed, and the simulation resolution accuracy is reduced.

The traditional random natural fracture modeling is based on near-wellbore geostatistics and well logging parameter fracture interpretation (tendency, inclination angle, trend and the like), local natural fracture description distortion is easily caused by adopting interpolation or random simulation processing among wells, particularly, only depending on a near-wellbore fracture density curve, realizing the random distribution processing hypothesis of natural fractures by using Monte Carlo and other statistical test methods, and coarsening the hypothesis to the whole grid to complete the modeling of the whole discrete fracture network, wherein the working process has the problems of low accuracy of prediction results, high uncertainty and the like. Secondly, for the simulation of the ground stress field of the special geological structure including faults, natural cracks and the like, the influences of distribution characteristics such as the faults and the cracks and mechanical properties of the distribution characteristics on the surrounding stress field need to be considered, and the influence of local stress characteristics on the global stress field is usually ignored in the existing ground stress modeling method, so that the final simulation result in the enriched region of the faults and the natural cracks is distorted.

The finite element method is adopted to simulate the distribution of the ground stress field, and the calculation result of the logging stress of a single well is usually used as a boundary condition, and an optimization algorithm is selected to realize iterative solution by using the boundary condition as a condition. However, the conventional particle swarm algorithm has limited initial samples, high probability of entering a local optimal solution, and difficulty in meeting the inverse regression of three-dimensional ground stress field data in efficiency and precision. A competitive particle swarm optimization (CSO) algorithm is one of evolutionary algorithms, and compared with a traditional algorithm, the CSO algorithm is simple in structure, good in performance and capable of being effectively applied to the optimal solution problem in an engineering scheme. However, the algorithm is not introduced into the numerical simulation of the three-dimensional ground stress field.

Competing particle populations differ from particle population optimization algorithms in that the optimal position and global optimal position of each particle in a competing particle population is no longer updated by pbest and gbest, but rather by a competing mechanism. The algorithm assumes that the size of a population is M, randomly initializes the population in a solution space, firstly randomly equally divides the population into 2 groups in each iteration process, two groups of particles compete and compare pairwise, and updates the information of a loser according to the information of a winner after each competition. The loser learns and updates his location and speed to the winner according to the following formula:

x_L(t+1)＝x_L(/)+v_L(t+1) (2)

in the formula: x is the number of_w(t)、x_L(t) position vectors representing winners and losers, respectively; v. of_L(t +1) represents the velocity vector of the loser; t is the number of iterations; r is₁(t)、r₂(t)、r₃(t)∈[0，1]^DThe method is characterized in that 3 random vectors which are subjected to uniform distribution have the same dimension with a solution vector; φ is a parameter that controls the effect of x (t) on loser location updates; x (t) has two meanings, one of which is the average position of all particles and is global; the other indicates the local average position of the particle in the domain, and has locality. The competitive particle swarm optimization solves the oscillation problem existing in the evolution of the conventional particle swarm optimization by mixing randomly selected targets in the self historical optimization and removing the global optimization, and improves the convergence precision and effect.

Through the above analysis, the problems and defects of the prior art are as follows:

(1) the prior art does not fully utilize seismic attribute information to constrain the prediction result among wells, particularly has limited resolution ratio in the aspect of fracture random modeling, and simultaneously does not analyze the influence of natural fracture development characteristics with different scales on a rock mechanical field and a ground stress field. In the traditional ground stress field simulation, the commonly adopted vertical main stress hypothesis and equal strain hypothesis sometimes cause that the self-consistency cannot be realized, and the simulation calculation precision of the ground stress field constructed by complex structures is low; for this reason, it is necessary to introduce a hierarchical load superposition and boundary cell extension technique to reduce the boundary effect of the ground stress model cell.

(3) In the prior art, the interaction between the geological corner grid and the finite element grid is realized by adopting a method of forced conversion or removing and sacrificing local node data, but the data is easy to lose and is abnormal due to the difference of the properties of the two grids. Therefore, node and unit information is searched by adopting a global and local double positioning and tracking method, node assignment is carried out by adopting a spherical region searching method and a natural adjacent point interpolation algorithm, the geological model attribute result is conducted to a stress model of a finite element, meanwhile, the method can also carry out inverse transformation of grid data, and important help is provided for the application of a subsequent ground stress field and the numerical reservoir simulation.

(3) In the prior art, a competitive particle swarm optimization algorithm is not introduced to solve and calculate the three-dimensional ground stress field, so that the obtained three-dimensional ground stress result is very low in precision and resolution. The convergence rate and the global search capability are two key criteria for measuring the performance of the particle swarm algorithm, and the competitive particle swarm algorithm can balance the convergence rate and the global search capability of particles and effectively improve convergence to global optimum.

The difficulty in solving the above problems and defects is: some technical difficulties exist in the process of processing the problems, such as how to extract relevant attributes through a seismic technology, perform deterministic modeling on natural fractures, and construct a rock mechanical model coupling the properties of the natural fractures on the basis, particularly, mutual influence among the natural fractures needs to be fully considered in the subsequent ground stress calculation, so that the distortion of the final stress field result obtained by a linear superposition method is avoided. Secondly, how to effectively position the nodes during the interaction of the grid data, and on the basis, the conversion rate and the precision between the grid attributes are influenced by searching the data according to the grid density and the node positions. In addition, when three-dimensional finite element iterative inversion is carried out, how to effectively control the deformation of a reservoir in a model space to ensure that the displacement and stress concentration of a grid unit gradually become smaller along with the distance from a research area and the displacement is zero and the original field stress level is reached on a boundary; secondly, the problems of result oscillation and limited inversion precision exist in the conventional three-dimensional finite element numerical inversion process.

The significance of solving the problems and the defects is as follows: the problems that the trial calculation of the boundary load by adopting a general optimization method is slow in iterative convergence and limited in calculation precision can be solved. Meanwhile, the constructed three-dimensional ground stress field can provide important technical support for subsequent oil and gas field exploration and development projects such as dessert prediction, well structure design of drilling, fracturing stimulation scheme analysis and the like. The grid interaction scheme developed in the process of the ground stress modeling is beneficial to the joint use of a plurality of simulators such as geology-oil deposit-fracturing and the like in the follow-up process, and the engineering application value of the simulators is improved.

Disclosure of Invention

In order to overcome the problems in the related art, the disclosed embodiment of the invention provides a reservoir three-dimensional stress field simulation method, a reservoir three-dimensional stress field simulation system, a reservoir three-dimensional stress field simulation terminal and a reservoir three-dimensional stress field simulation storage medium. In particular to a reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization.

The technical scheme is as follows: according to the first aspect provided by the invention, a reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization is provided, and comprises the following steps:

classifying the seismic facies according to the development characteristics of bedding cracks and natural cracks on the well;

under the constraint of the seismic facies, the crack lines are tracked by using the enhanced curvature attribute, and a deterministic model of a discrete natural crack network is established on the seismic scale;

based on the established fracture certainty model, aiming at the influence of local weakening such as faults, natural fractures and the like on the whole rock mechanical field and the ground stress field, the grid conversion of the geological model and the ground stress finite element model is carried out by utilizing a related grid conversion program, and the optimal solution acquisition of the ground stress field is realized by utilizing a competitive particle swarm optimization algorithm.

In an embodiment of the present invention, the reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization specifically includes:

step one, building a block structure model by using a well-seismic combination mode, comprising the following steps of: forming a whole-region construction model by using a two-dimensional seismic interpretation result as a transverse basis, using single-well layered data as a correction point and combining fault development characteristics;

determining longitudinal and transverse constraint conditions and selecting a variation function, and performing lithofacies modeling by using multi-point geostatistics and a sequential indication random simulation method; analyzing distribution ranges and frequency distribution forms of various attributes of different lithofacies on the basis of lithofacies modeling and structural modeling, performing key attribute modeling such as permeability, porosity and the like by an interpolation method, and performing constraint by single well interpretation data;

step three, utilizing a three-dimensional earthquake to describe and depict natural cracks of the block;

preprocessing a logging curve, correcting logging data of the expanding and well wall scouring positions, wherein the logging data comprise stratum gamma rays, density and longitudinal wave time, and rock parameters obtained through calculation according to acoustic logging data comprise elastic modulus and Poisson ratio key parameters; performing a uniaxial/triaxial compression test and a longitudinal and transverse sound wave speed test on the mechanical parameters of the rock core of the target layer, and analyzing and establishing a dynamic and static parameter conversion model;

step five, obtaining dynamic elastomer parameters through full-time window prestack three-dimensional seismic inversion from the earth surface to a target layer, calculating static data by using a dynamic and static conversion model to obtain a three-dimensional rock mechanics data body, and introducing geological statistical information to predict thin layer elasticity parameters for the thin layer mechanics parameters;

step six, introducing the fault serving as a discontinuous body into a three-dimensional rock mechanical parameter model, and analyzing the variation of the ground stress field at two sides of the fault; correcting the three-dimensional rock mechanical parameter model by combining single well logging, seismic interpretation results and indoor experimental results, and finally updating the three-dimensional rock mechanical field combined with well and seismic to obtain a three-dimensional rock mechanical parameter model aiming at geological structure influence;

step seven, utilizing single well acoustic-resistivity and density logging data to predict pore pressure, utilizing on-site drilling or logging actual measurement data to constrain, integrating drilling and well testing data, determining the interval between the maximum and minimum pore pressure values of a target layer, and constraining the pore pressure calculation result;

step eight, combining the characteristics of block geological stratification, reversely solving by using field collected small pressure data and well drilling data to obtain the structural stress coefficient of each layer, carrying out single-well logging ground stress calculation by combining a spring model ground stress calculation formula or other ground stress formulas, and correcting the logging calculation ground stress data by using the ground stress experiment results (difference strain experiment or Kessel acoustic emission experiment) of different layers;

step nine, converting the geological corner point grid into a corresponding three-dimensional finite element grid model by using a grid interaction algorithm, searching node and unit information by using a global and local positioning method, and performing node assignment by using a spherical region searching method and a natural neighbor interpolation algorithm to realize the transmission of the geological model attribute result to a stress model of a finite element;

step ten, dispersing the block geometry;

step eleven, after determining a ground stress inversion work area and parameters of participation, establishing a three-dimensional ground stress field inversion target function, initializing particle swarm parameters, determining errors and iteration precision, and obtaining a global optimal solution by utilizing a competitive particle swarm algorithm;

and step twelve, carrying out finite element numerical simulation calculation, stopping calculation after iteration or calculation precision, and finally obtaining the three-dimensional ground stress distribution characteristics of the block.

In an embodiment of the present invention, the third step of describing and describing the block natural fracture by using a three-dimensional earthquake specifically includes:

acquiring distribution characteristics of natural fractures by collecting multi-source information obtained by field outcrop, single-well core observation, logging cuttings and early geological research results in a research area, and constructing a bedding crack model by taking the established rock phase model and thickness reservoir model as constraints by combining with reservoir bedding weak plane explanation;

carrying out deterministic modeling on dispersion cracks by adopting a grouping mode to quantify the crack properties of the small-scale dispersion cracks;

performing topological parameters by adopting an ant colony tracking algorithm, establishing a related natural fracture network for large-scale discrete fracture analysis, forming a three-dimensional seismic data body by a plurality of two-dimensional seismic data surfaces, performing grid division on the two-dimensional seismic data surfaces and defining ant initial probability;

wherein, P_iIs a probability, C_iIs the coherence value at the ith point; under the condition of determining the starting probability, determining a main tracking direction by further dividing the region and calculating the gradient; the trace direction of each sub-block is expressed as

Wherein, theta_iTo track the direction G_x(i, j) and G_y(i, j) is the gradient of the data points.

In an embodiment of the present invention, the linear formula of the prestack elastic parameter inversion of the reflectivity of the step five based on direct inversion of young's modulus, poisson's ratio and density is as follows:

in the formula: r is a reflection coefficient; theta is an incident angle; v. of_PIs the velocity of the longitudinal wave; v. of_sIs the transverse wave velocity; k is

E is Young's modulus; sigma is the Poisson ratio; ρ is the density.

In an embodiment of the invention, in the sixth step, a mechanical parameter field after coupling the natural fracture is established for the influence of the fracture stiffness on the rock mechanical parameters;

Δμ_n: normal displacement, mm, k_n: normal stiffness, MPa/mm, k_s: shear stiffness, MPa/mm

Δτ_s: change in shear stress, MPa; Δ μ_s: tangential displacement, mm; delta tau_n: change in normal stress, MPa

The normal stress of the crack is in direct proportion to the rigidity of the crack, and the size of the crack is in inverse proportion to the rigidity of the crack; the rigidity parameter of the crack is drawn up according to the Bandis experimental data;

fracture shear stiffness is a function related to scale, fracture face normal stress, as follows:

in the formula: k_sShear stiffness of crack, σ_nThe normal stress of the crack surface, r the crack size, JRC the roughness, JCS the compressive strength,

is the rubbing angle.

In an embodiment of the invention, the seventh step is to use the measured data of the on-site drilling or logging to carry out constraint, and if the description of blowout, kick and obvious single gas does not appear in the drilling process, the drilling mud density is used as the upper limit of the formation pressure; and (3) balancing the stratum static pressure obtained in the well test with the stratum fluid in the shaft to obtain the stratum pressure, and if the balancing time is insufficient, taking the obtained static pressure as the lower limit of the stratum pressure.

In an embodiment of the present invention, in the step eight, a calculation formula of the combined spring model horizontal stress model is as follows:

in the ninth step, a finite element load method is used for calculating a ground stress field, the three-dimensional research area is expanded in the horizontal direction, boundary grids are added, and the geometric size of each side of the expanded grid is 2-3 times that of the research area; meanwhile, a virtual grid is added, the three-dimensional grid extends upwards to the ground surface, and extends downwards to a proper depth below the reservoir stratum, so that the stress strain in the research area meets the continuous condition; meanwhile, introducing a new layer in the lowest layer of the stress model of the finite element, setting the bottom surface as a plane, and adding displacement constraint; taking the structure, the large fault and the fracture as stress model geometric boundaries of finite elements, and combining the gravity action to obtain the ground stress distribution of a horizontal stress boundary and a horizontal displacement boundary under the action of a load;

the step ten specifically comprises: the finite element grid adopts a random search method to carry out quality correction and treatment on the distorted and elongated grid so that the physical and mechanical properties of each rock are smooth and continuous between the grid and the grid, and finally, the amplitude of the parameter of the stratum material is obtained through well logging calculation data and well drilling actual measurement information.

According to a second aspect provided by the invention, a reservoir three-dimensional stress field simulation system based on deterministic multi-scale fracture modeling and competitive particle swarm optimization comprises:

the three-dimensional geological modeling module is used for obtaining related geological structures by interpreting seismic wave data of the whole bed system, correcting the seismic data by using single well data and geological stratification data, establishing an accurate block structure model by combining well-seismic information, modeling geological attributes by using a well-seismic combination modeling technology on the basis of three-dimensional structural modeling and phase modeling, and dividing fine layers of regions and platforms;

the natural fracture modeling module is used for constructing a bedding fracture model by taking the established rock phase model and the thickness reservoir model as constraints and combining with reservoir bedding weak plane explanation based on field outcrop, single-well core analysis, logging rock debris and early-stage geological analysis results of a research area; constructing a natural fracture model by combining the occurrence, density and development influence factors of the structural fractures in the research area, tracking fracture lines by using enhanced curvature attributes under the constraint of the seismic facies, realizing the establishment of a deterministic model of a discrete natural fracture network on the seismic scale, and correcting and debugging the natural fracture modeling result by combining near-wellbore geostatistics and well logging parameter fracture interpretation;

the logging curve correction and quality analysis module is used for performing single-well rock mechanical parameter and ground stress explanation by using a rock mechanical and ground stress model formula, and performing single-well geomechanical calculation result restriction by using relevant data such as a single-shaft/three-shaft compression test, a longitudinal and transverse sound wave speed test and the like of the rock mechanical parameter development;

the three-dimensional attribute model module is used for establishing a longitudinal and transverse wave and density three-dimensional attribute model of a research area by utilizing pre-stack seismic inversion, calculating the Young modulus and Poisson ratio of the rock, and obtaining a three-dimensional rock mechanical model on the basis of the distribution trend of seismic attribute bodies;

the three-dimensional rock mechanical field model correction module of the natural fracture is used for assigning measurement parameter values and model parameter initial values to different layers of the mechanical model according to the modeling results of the natural fracture and the fault and by combining the layer information of a single well; taking the fault as a discontinuous body, introducing the fault into a three-dimensional rock mechanical field model, and analyzing the change of the ground stress field at two sides of the fault; the method is also used for calculating the influence of the fracture rigidity on the rock mechanics parameter and establishing a three-dimensional rock mechanics parameter field after coupling the natural fracture;

the single-well stratum pore pressure model module is used for analyzing the deposition rule and the abnormal pressure generation mechanism of the research area, comprehensively utilizing logging data such as sound wave-resistivity and density to establish a single-well stratum pore pressure model, utilizing well drilling and completion and oil testing data to carry out data correction and quality control to obtain the relationship between the stratum pressure coefficient and the burial depth, lithology and porosity, and establishing the single-well stratum pore pressure model of the research area;

the grid conversion algorithm application module is used for searching node and unit information by adopting a global and local double positioning and tracking method, performing node assignment by adopting a spherical region searching method and a natural neighbor interpolation algorithm, converting geological corner grids into corresponding three-dimensional finite element grid models, and simultaneously utilizing a grid conversion reverse method and resampling the models to the corresponding corner grid geological models;

the rock stratum boundary unit adding module is used for adding overlying rock stratum, underlying rock stratum and lateral rock stratum boundary units for the research area, controlling the deformation of the reservoir stratum in a model space, and gradually reducing the displacement and stress concentration of the grid units along with the distance from the research area to reach the level of zero displacement and original site stress on the boundary;

the local ground stress field direction acquisition module is used for collecting drilling induced joints and well bore caving fracture trends of different layers of a work area, and forming an imaging rose diagram of each layer by combining FMI imaging logging results to obtain a local ground stress field direction;

the inversion grid boundary module is used for setting corresponding load and constraint conditions for a target reservoir, adding attributes for related media of the finite element model, performing elastoplasticity finite element calculation, and inverting the constraint, stress magnitude and direction conditions of the grid boundary;

and the competitive particle swarm algorithm solving module is used for acquiring a global optimal solution and a particle historical optimal solution by using the competitive particle swarm algorithm, comparing a calculation result with an actual measurement result of the ground stress size and direction of a known point, and stopping calculation after the accuracy or the iteration times are met.

According to a third aspect provided by the present invention, there is provided an information data processing terminal comprising a memory and a processor, the memory storing a computer program which, when executed by the processor, causes the processor to execute the reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization.

According to a fourth aspect provided by the present invention, there is provided a computer readable storage medium storing a computer program which, when executed by a processor, causes the processor to execute the reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization.

By combining all the technical schemes, the invention has the advantages and positive effects that:

the invention provides a prediction method for carrying out deterministic natural fracture modeling and ground stress modeling by combining well earthquake and performing numerical simulation inversion of a ground stress field by using a competitive particle swarm algorithm on the basis, which is mainly innovative in classifying earthquake phases according to development characteristics of bedding fractures and natural fractures on a well and hydraulic fracturing degree. Under the constraint of the seismic facies, the crack lines are tracked by utilizing the enhanced curvature property, a deterministic model of a discrete natural crack network is established on the seismic scale, the influence of local weakening such as faults, natural cracks and the like on the whole rock mechanical field and the ground stress field is considered, a related grid conversion program is developed to realize grid conversion of a geological model and a ground stress finite element model, and the optimal solution acquisition of the ground stress field is realized by utilizing a competitive particle swarm optimization algorithm.

The invention improves the accuracy of natural fracture modeling, considers the influence of natural fracture development characteristics of different scales on a rock mechanical field and a ground stress field, effectively solves the problems that vertical principal stress hypothesis, equal strain hypothesis and the like cannot be consistent, utilizes seismic information to constrain inter-well prediction results, adopts global and local double positioning and tracking methods to search node and unit information, and adopts a spherical region search method and a natural neighbor interpolation algorithm to carry out node assignment, realizes the conversion and interaction of multi-model unit grids and attributes, proposes that in rock mechanical modeling, the influence of local strength weakening caused by geological discontinuous surfaces such as faults and natural fractures on the rock mechanical field is considered, ground stress inversion is carried out on the basis, simultaneously a competitive particle swarm optimization algorithm is introduced to realize three-dimensional ground stress field solving calculation, compared with the existing linear superposition ground stress inversion technology, the simulation calculation precision and speed of the complex structure ground stress field can be greatly improved, and the three-dimensional ground stress result with high precision and resolution can be finally obtained, so that important technical support is provided for subsequent oil and gas field exploration and development work such as engineering dessert prediction, well structure design of a drilling well, fracturing yield-increasing scheme analysis and the like.

It is to be understood that both the foregoing general description and the following detailed description are exemplary and explanatory only and are not restrictive of the disclosure.

Drawings

The accompanying drawings, which are incorporated in and constitute a part of this specification, illustrate embodiments consistent with the present disclosure and together with the description, serve to explain the principles of the disclosure.

Fig. 1 is a flow chart of a reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization provided by the embodiment of the invention.

Fig. 2 is a schematic diagram of a three-dimensional ground stress field modeling method provided by an embodiment of the invention.

FIG. 3 is a graph of data from the application of natural fracture modeling results provided by embodiments of the present invention.

Fig. 4 is a graph of application data of rock mechanics numerical simulation considering natural fracture weakening provided by an embodiment of the present invention. FIG. 4 (a); FIG. 4 (b).

Fig. 5 is a graph of the result of the three-dimensional ground stress field modeling application provided by the embodiment of the invention. FIG. 5(a) minimum horizontal principal stress distribution plot; (b) maximum horizontal principal stress profile.

FIG. 6 is a graph comparing the three-dimensional ground stress results and the one-dimensional ground stress magnitude results provided by the embodiment of the present invention.

Detailed Description

In order to make the aforementioned objects, features and advantages of the present invention comprehensible, embodiments accompanied with figures are described in detail below. In the following description, numerous specific details are set forth in order to provide a thorough understanding of the present invention. This invention may, however, be embodied in many different forms and should not be construed as limited to the embodiments set forth herein, but rather should be construed as broadly as the present invention is capable of modification in various respects, all without departing from the spirit and scope of the present invention.

It will be understood that when an element is referred to as being "secured to" another element, it can be directly on the other element or intervening elements may also be present. When an element is referred to as being "connected" to another element, it can be directly connected to the other element or intervening elements may also be present. As used herein, the terms "vertical," "horizontal," "left," "right," and the like are for purposes of illustration only and are not intended to represent the only embodiments.

Unless defined otherwise, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention belongs. The terminology used in the description of the invention herein is for the purpose of describing particular embodiments only and is not intended to be limiting of the invention. As used herein, the term "and/or" includes any and all combinations of one or more of the associated listed items.

As shown in FIG. 1, the invention provides a reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization, which comprises the following steps:

s101, three-dimensional geological modeling, namely firstly, obtaining related geological structures by using the explained seismic wave data of the whole bed system, correcting the seismic data by using single well data and geological stratification data, establishing an accurate block structure model by combining well-seismic information, then, performing geological attribute modeling by using a well-seismic combination modeling technology on the basis of three-dimensional structural modeling and phase modeling, and carrying out fine sub-layer division of regions and platforms.

S102, natural fracture modeling, namely constructing a bedding fracture model by taking the established rock facies model and the thickness reservoir model as constraints and combining with reservoir bedding weak plane explanation based on field outcrop of a research area, single-well core observation, logging rock debris and early geological research results. And (2) constructing a natural fracture model by combining the occurrence and density of the structural fractures in the research area and development influence factors (such as distance from a fault and the like), tracking fracture lines by utilizing an enhanced curvature attribute under the constraint of the seismic facies, establishing a deterministic model of a discrete natural fracture network on the seismic scale, and correcting and debugging the modeling result of the natural fractures by combining near-wellbore geostatistics and well logging parameter fracture interpretation (tendency, inclination angle, trend and the like).

S103, well logging curve correction and quality analysis are carried out, single-well rock mechanical parameters and ground stress interpretation is carried out by utilizing a rock mechanical and ground stress model formula, and single-well geomechanical calculation results are restrained by utilizing relevant data such as a single-shaft/three-shaft compression test, a longitudinal and transverse sound wave speed test and the like of rock core mechanical parameters.

S104, establishing a three-dimensional attribute model of longitudinal and transverse waves, density and the like in a research area by utilizing pre-stack seismic inversion, calculating the Young modulus and Poisson ratio of the rock, and obtaining a three-dimensional rock mechanical model on the basis of the distribution trend of seismic attribute bodies.

And S105, combining the correction of the three-dimensional rock mechanical field model of the natural fracture, and assigning measurement parameter values and model parameter initial values to different layers of the mechanical model according to the modeling results of the natural fracture and the fault and combining the layer position information of a single well. The fault is considered as a discontinuous body and introduced into a model, and the change of the ground stress field at two sides of the fault is analyzed. The analysis of the fracture zone is similar to that of the fault, and the fracture zone is affected more by rock mechanics parameters than the developed enrichment zone. Secondly, calculating the influence of the fracture rigidity on the rock mechanics parameter by using a correlation formula, and establishing a three-dimensional rock mechanics parameter field coupled with the natural fracture.

S106, analyzing the deposition law and the abnormal pressure generation mechanism of the research area, comprehensively utilizing logging data such as sound wave-resistivity and density to establish a single-well stratum pore pressure model, utilizing well drilling and completion and oil testing data to carry out data correction and quality control, obtaining the relation between a stratum pressure coefficient and the burial depth, lithology and porosity, and establishing the single-well stratum pore pressure model of the research area.

S107, a grid conversion algorithm is applied, node and unit information is searched by adopting a global and local double positioning and tracking method, node assignment is carried out by adopting a spherical region searching method and a natural neighbor interpolation algorithm, geological corner grids are converted into corresponding three-dimensional finite element grid models, and meanwhile, the models can be resampled to the corresponding corner grid geological models by utilizing a grid conversion reverse method.

And S108, adding boundary units of an overlying rock stratum, a lower overlying rock stratum and a lateral rock stratum into the research area, wherein the purpose is mainly to control the deformation of the reservoir in a model space, and in addition, the displacement and stress concentration of the grid units are gradually reduced along with the distance from the research area, so that the displacement and the original field stress level are reached on the boundary.

S109, collecting drilling induced joints and well bore caving fracture trends of different layers of a work area, and forming an imaging rose diagram of each layer by combining FMI imaging logging results to obtain the local ground stress field direction.

S110, setting corresponding load and constraint conditions for a target reservoir, adding attributes for relevant media of a finite element model, performing elastoplasticity finite element calculation, inverting the constraint and stress magnitude and direction conditions of a grid boundary, simultaneously obtaining a global optimal solution and a particle historical optimal solution by using a competitive particle swarm algorithm, comparing a calculation result with the ground stress magnitude and direction actual measurement result of a known point, and stopping calculation after the precision or iteration times are met.

The invention also provides a reservoir three-dimensional stress field simulation system based on deterministic multi-scale fracture modeling and competitive particle swarm optimization, which comprises the following steps:

the three-dimensional geological modeling module is used for obtaining related geological structures by interpreting seismic wave data of the whole bed system, correcting the seismic data by using single well data and geological stratification data, establishing an accurate block structure model by combining well-seismic information, modeling geological attributes by using a well-seismic combination modeling technology on the basis of three-dimensional structural modeling and phase modeling, and dividing fine layers of regions and platforms;

the natural fracture modeling module is used for constructing a bedding fracture model by taking the established rock phase model and the thickness reservoir model as constraints and combining with reservoir bedding weak plane explanation based on field outcrop, single-well core analysis, logging rock debris and early-stage geological analysis results of a research area; constructing a natural fracture model by combining the occurrence, density and development influence factors of the structural fractures in the research area, tracking fracture lines by using enhanced curvature attributes under the constraint of the seismic facies, realizing the establishment of a deterministic model of a discrete natural fracture network on the seismic scale, and correcting and debugging the natural fracture modeling result by combining near-wellbore geostatistics and well logging parameter fracture interpretation;

the logging curve correction and quality analysis module is used for performing single-well rock mechanical parameter and ground stress explanation by using a rock mechanical and ground stress model formula, and performing single-well geomechanical calculation result restriction by using relevant data such as a single-shaft/three-shaft compression test, a longitudinal and transverse sound wave speed test and the like of the rock mechanical parameter development;

the three-dimensional attribute model module is used for establishing a longitudinal and transverse wave and density three-dimensional attribute model of a research area by utilizing pre-stack seismic inversion, calculating the Young modulus and Poisson ratio of the rock, and obtaining a three-dimensional rock mechanical model on the basis of the distribution trend of seismic attribute bodies;

the three-dimensional rock mechanical field model correction module of the natural fracture is used for assigning measurement parameter values and model parameter initial values to different layers of the mechanical model according to the modeling results of the natural fracture and the fault and by combining the layer information of a single well; taking the fault as a discontinuous body, introducing the fault into a three-dimensional rock mechanical field model, and analyzing the change of the ground stress field at two sides of the fault; the method is also used for calculating the influence of the fracture rigidity on the rock mechanics parameter and establishing a three-dimensional rock mechanics parameter field after coupling the natural fracture;

the single-well stratum pore pressure model module is used for analyzing the deposition rule and the abnormal pressure generation mechanism of the research area, comprehensively utilizing logging data such as sound wave-resistivity and density to establish a single-well stratum pore pressure model, utilizing well drilling and completion and oil testing data to carry out data correction and quality control to obtain the relationship between the stratum pressure coefficient and the burial depth, lithology and porosity, and establishing the single-well stratum pore pressure model of the research area;

the grid conversion algorithm application module is used for searching node and unit information by adopting a global and local double positioning and tracking method, performing node assignment by adopting a spherical region searching method and a natural neighbor interpolation algorithm, converting geological corner grids into corresponding three-dimensional finite element grid models, and simultaneously utilizing a grid conversion reverse method and resampling the models to the corresponding corner grid geological models;

the rock stratum boundary unit adding module is used for adding overlying rock stratum, underlying rock stratum and lateral rock stratum boundary units for the research area, controlling the deformation of the reservoir stratum in a model space, and gradually reducing the displacement and stress concentration of the grid units along with the distance from the research area to reach the level of zero displacement and original site stress on the boundary;

the local ground stress field direction acquisition module is used for collecting drilling induced joints and well bore caving fracture trends of different layers of a work area, and forming an imaging rose diagram of each layer by combining FMI imaging logging results to obtain a local ground stress field direction;

the inversion grid boundary module is used for setting corresponding load and constraint conditions for a target reservoir, adding attributes for related media of the finite element model, performing elastoplasticity finite element calculation, and inverting the constraint, stress magnitude and direction conditions of the grid boundary;

and the competitive particle swarm algorithm solving module is used for acquiring a global optimal solution and a particle historical optimal solution by using the competitive particle swarm algorithm, comparing a calculation result with an actual measurement result of the ground stress size and direction of a known point, and stopping calculation after the accuracy or the iteration times are met.

The technical solution of the present invention is further described below with reference to specific examples.

Example, as shown in fig. 2.

The invention provides a reservoir three-dimensional stress field simulation method based on deterministic modeling and deep learning. The method comprises the following steps:

step 1: a block structure model is established by utilizing a well-seismic combination mode, namely a two-dimensional seismic interpretation result is used as a transverse basis, single-well layered data is used as a correction point, and a fault development characteristic is combined to form a whole-area structure model.

Step 2: and (3) considering the longitudinal distribution rule of the lithofacies, determining longitudinal and transverse constraint conditions and selecting a variation function, and performing lithofacies modeling by using geostatistics data and a sequential indication random simulation method. On the basis of lithofacies modeling and structural modeling, the distribution range and the frequency distribution form of each attribute of different lithofacies are analyzed, key attributes such as permeability, porosity and the like are modeled through an interpolation method, and constraint is performed through single well interpretation data.

And step 3: and (4) utilizing a three-dimensional earthquake to depict and describe the natural cracks of the block. Among the deterministic inversion methods, constrained sparse pulse inversion is the most important method. The constraint sparse pulse inversion is based on inversion of the traces, and the essence of the inversion is that the optimal matching of the synthetic record and the seismic traces is achieved by using the minimum number of reflection coefficient pulses under the constraint of the impedance trend. The whole natural fracture is divided into three types, one type is a bedding fracture, the transverse isotropy of rock mechanics is usually determined by the fractures, the accurate acquisition of the distribution characteristics of the fractures is realized by collecting multi-source information obtained by field outcrop, single-well core observation, well logging cuttings and early geological research results in a research area, and meanwhile, the construction of a bedding fracture model is realized by taking the established rock facies model and the thickness reservoir model as constraints by combining with reservoir bedding weak face explanation. The second type is small-scale dispersion fracture, and in practical research, the maximum curvature and the azimuth angle attribute thereof are used as the development density and trend of the dispersion fracture. In order to better quantify the fracture properties, deterministic modeling of dispersion fractures is achieved in a grouping manner. The third type is large-scale discrete fracture, which determines the connectivity of the entire fracture. For the cracks, an ant colony tracking algorithm is mainly adopted to carry out topological parameters, and a related natural crack network is established. The three-dimensional seismic data volume is essentially composed of a plurality of two-dimensional seismic data surfaces, so that the two-dimensional seismic data surfaces can be subjected to grid division and the ant starting probability can be defined.

Wherein, P_iIs a probability, C_iIs the coherence value at the ith point. On the premise that the starting probability is determined, the main direction of tracking is determined by further dividing the region and calculating the gradient. Tracing direction of each sub-blockCan be expressed as

Wherein, theta_iTo track the direction G_x(i, j) and G_y(i, j) is the gradient of the data points.

And 4, step 4: and preprocessing a logging curve, correcting logging data of the expanding and well wall scouring positions, wherein the logging data comprise the gamma rays, the density and the longitudinal wave duration of the stratum, and rock parameters obtained by calculation according to the acoustic logging data comprise key parameters such as the elastic modulus, the Poisson ratio and the like. And (3) carrying out a uniaxial/triaxial compression test and a longitudinal and transverse sound wave speed test on the mechanical parameters of the rock core of the target layer, and analyzing and establishing a dynamic and static parameter conversion model.

And 5: considering the problems of low frequency, lack of structural information and the like of conventional seismic speed, the full-range system carries out high-precision inversion to finely depict elastic information, layered stratum speed details are completely presented, dynamic elastic body parameters are obtained through full-time window pre-stack three-dimensional seismic inversion from the earth surface to a target layer, static data calculation is carried out by utilizing a dynamic and static conversion model, a three-dimensional rock mechanics data body is obtained, and for the thin layer mechanics parameters, geological statistical information needs to be introduced to realize thin layer elastic parameter prediction. The linear formula of the prestack elastic parameter inversion of the reflectivity based on the direct inversion of the Young modulus, the Poisson ratio and the density is as follows:

in the formula: r is a reflection coefficient; theta is an incident angle; v. of_PIs the velocity of the longitudinal wave; v. of_sIs the transverse wave velocity; k is

E is Young's modulus; sigma is the Poisson ratio; ρ is the density.

Step 6: and (4) introducing the fault into the model as a discontinuous body, and analyzing the change of the ground stress field at two sides of the fault. The analysis of the fracture zone is similar to that of the fault, and the fracture zone is affected more by rock mechanics parameters than the developed enrichment zone. And (4) considering the influence of the fracture rigidity on the mechanical parameters of the rock mass, and establishing a mechanical parameter field after coupling the natural fracture.

Δμ_n: normal displacement, mm, k_n: normal stiffness, MPa/mm, k_s: shear stiffness, MPa/mm

Δτ_s: change in shear stress, MPa; Δ μ_s: tangential displacement, mm; delta tau_n: change in normal stress, MPa

The normal stress of the fracture is proportional to the stiffness of the fracture, and the size of the fracture is inversely proportional to the stiffness of the fracture. The stiffness parameters of the crack can be formulated according to experimental data of Bandis and the like. The fracture shear stiffness is a function related to the dimension and the normal stress of the fracture surface, and the basic principle of the Oda method is as follows:

in the formula: k_sShear stiffness of crack, σ_nThe normal stress of the crack surface, r the crack size, JRC the roughness, JCS the compressive strength,

is the rubbing angle.

And correcting the three-dimensional rock mechanical parameter model by combining single well logging, seismic interpretation results and indoor experimental results, and finally updating the three-dimensional rock mechanical field combined with well and seismic to obtain the three-dimensional rock mechanical parameter model considering the influence of geological structures.

And 7: the pore pressure prediction is carried out by using single-well acoustic-resistivity and density logging data, preferably using an Eton method, a Bowers effective stress method and the like as calculation methods, and using field drilling or logging and other measured data for constraint, wherein if well blowout, well kick and description of obvious single gas do not occur in the drilling process, the drilling mud density can be used as the upper limit of the formation pressure (coefficient). The static pressure of the stratum obtained in the well test can be used as the lower limit of the stratum pressure (especially in low-permeability stratum) if the balance between the borehole and the stratum fluid is achieved, and the stratum pressure can be obtained if the balance time is insufficient. And by integrating the drilling and well testing data, the maximum and minimum range intervals of the pore pressure value of the target layer can be determined, so that the calculation result of the pore pressure is restrained.

And 8: meanwhile, combining the characteristics of block geological stratification, reversely solving by using field collected small pressure data and well drilling data to obtain the structural stress coefficient of each layer, calculating the ground stress of single-well logging by using the following combined spring model formula, and correcting the logging data by using the ground stress calculation results of different layers. The combined spring model horizontal stress model has the calculation formula as follows:

and step 9: and converting the geological corner point grid into a corresponding three-dimensional finite element grid model by utilizing a grid interaction algorithm, and transmitting the geological model attribute result into a stress model of a finite element. And (3) carrying out ground stress field calculation by using a finite element load method, expanding the three-dimensional research area in the horizontal direction, adding boundary grids, wherein the geometric dimension of each expanded grid is 2-3 times of that of the research area. Meanwhile, a virtual grid is properly added, the three-dimensional grid extends upwards to the ground surface, and the downward edge extends to a proper depth below the reservoir, so that stress and strain in the research area are ensured to meet continuous conditions, and the whole model is consistent and stable. Meanwhile, a new horizon is introduced into the lowest layer of the model, and the bottom surface is set as a plane, so that displacement constraint is added. And comprehensively considering the gravity action by taking the structure, the large fault and the fracture as the geometric boundary of the model body, and acquiring the ground stress distribution of the model under the estimation load action by the horizontal stress boundary and the horizontal displacement boundary.

Step 10: the block geometry is discretized. On the premise of considering the resolution of input/output data, the number of the three-dimensional geomechanical finite element grid systems is determined to ensure the calculation precision and speed, and the grid scale is controlled to be in the level of ten million. The finite element grid should be as regular as possible, and the distorted and elongated grid is subjected to quality correction and processing by adopting a random search method, so that the transition of the physical and mechanical properties of each rock between the grid and the grid is ensured to be more smooth and continuous, and the deformation characteristics and the solving precision of the reservoir are met. And finally, the amplitude of the formation material parameter is measured through the well logging calculation data and the well drilling actual measurement information.

Step 11: after determining a ground stress inversion work area and parameters of participation, establishing a three-dimensional ground stress field inversion target function, initializing particle swarm parameters, determining errors and iteration precision, and obtaining a global optimal solution by utilizing a competitive particle swarm algorithm.

Step 12: and carrying out finite element numerical simulation calculation, stopping calculation after iteration or calculation precision, and finally obtaining the three-dimensional ground stress distribution characteristics of the block.

The technical solution of the present invention will be further described with reference to the positive effects.

The invention is based on the knowledge of sedimentary laws, and utilizes seismic structure interpretation and layered data and a well-seismic combination mode to establish a block structure model, namely, a two-dimensional seismic interpretation result is used as a transverse basis, single-well layered data is used as a correction point, and fault development characteristics are combined to form a whole-region structure framework.

The method uses interwell reservoir contrast to count the reservoir extension length as the main basis of the variation range of the model variation function, uses the reservoir thickness distribution diagram as the model plane trend constraint, uses the lithofacies vertical probability distribution as the model vertical trend constraint, and establishes the lithofacies model containing the reservoir layer through a sequential indication simulation method. The horizontal trend surface constraint of the lithofacies model mainly adopts a reservoir sand body thickness map made in a geological research stage, and the lithological character of a few-well or no-well area can be well controlled to be predicted, so that the horizontal consistency of model prediction and geological rules is realized. And (3) explaining the spatial correlation of reservoir parameters by using a variation function, and realizing lithofacies modeling through the longitudinal and transverse constraint conditions.

The method is used for counting the difference of different lithologies on attributes, the distribution range of each attribute and the frequency distribution form, then performing oil reservoir attribute modeling by using a well-seismic combined modeling technology on the basis of three-dimensional structural modeling and phase modeling, and performing multiple regression acquisition on key attributes such as gas-containing or oil-containing saturation by using measured data. The porosity and permeability properties need to be combined with field midway tests or petrophysical experiments to obtain the regular distribution of the porosity and permeability properties.

The whole natural fracture network is divided into a bedding fracture, a small-scale dispersion fracture, a large-scale dispersion fracture network and a fracture corridor. The bedding cracks and small-scale dispersed cracks are generally distributed in the whole oil reservoir, the effect of remarkably enhancing the communication capacity among the dispersed cracks is achieved in the fracturing modification, and the distribution of the bedding cracks and the small-scale dispersed cracks in the oil and gas reservoir is controlled by the mechanical property, the physical property and the structural position of rocks. For bedding joint modeling, the established rock facies model and thickness reservoir model are used as constraints in combination with multi-source information obtained by field outcrop of a research area, single-well core observation, logging rock debris and early geological research results, and meanwhile reservoir bedding weak plane explanation is utilized. For small-scale dispersion and large-scale dispersion cracks, the matching synergistic relation of time and space needs to be considered, seismic geometric multiple attributes under different resolution scales are proposed, and seismic facies are classified according to the development characteristics and the fracturing degree of aboveground cracks. Under the constraint of seismic facies, the crack lines are tracked by utilizing the enhanced curvature attribute, and the deterministic model of the discrete natural crack network is established on the seismic scale.

The method utilizes the acoustic logging data to calculate the mechanical parameters of the single-well rock, including the key parameters such as elastic modulus, Poisson ratio and the like. And (3) carrying out uniaxial/triaxial compression test and longitudinal and transverse sound wave speed test on the mechanical parameters of the rock core of the target layer, and simultaneously establishing a dynamic and static parameter conversion model for analysis.

According to the invention, complete elastic parameters are obtained through a full-layer prestack inversion technology, and the complete rock elastic parameters based on the seismic velocity are obtained through calculation.

The method considers the influence of the fracture rigidity on the rock mechanical parameters, assigns the mechanical parameters of the natural fracture and the fault position, and establishes the mechanical parameter field after coupling the natural fracture.

According to the method, a three-dimensional rock mechanical parameter model is corrected by combining single-well logging, seismic interpretation results and indoor experimental results, and finally, a well-seismic combined three-dimensional rock mechanical field is updated to obtain the three-dimensional rock mechanical parameter model considering geological structure influence.

The geological model attribute result is transmitted to a stress model of a finite element, and boundary units of an overlying rock stratum, an underlying rock stratum and a lateral rock stratum are added, wherein a core reservoir layer, the upper stratum and the lower stratum of the core reservoir layer are dispersed by adopting a 20-node second-order unit, and the rest of the strata are dispersed by adopting an 8-node linear unit. And adjusting the vertical height of the grid, and carrying out grid encryption on the core layer position by using a local encryption method so as to conveniently obtain a calculation result with higher resolution ratio in the subsequent process.

Load and material properties are added to the model, and the force borne by the model comprises three types: the method comprises the following steps of firstly, gravity load, secondly, pore pressure load and thirdly, boundary load. The boundary load is only loaded on the boundary of the model, and the direction and the magnitude of the stress in the model can change along with the local load and the rock mechanical property.

The method utilizes a competitive particle swarm algorithm to carry out ground stress inversion, firstly provides task set data, the maximum iteration times of the algorithm, a population rule and a variation step length, then carries out particle swarm initial chaotic initialization, calculates the particle fitness, carries out pairwise competition marking calculation on the whole population through the competitive particle swarm algorithm, sets a viewing command in the algorithm to carry out real-time acquisition on the algorithm completion, and can carry out parameter adjustment in real time when a program is premature or the algorithm is not convergent. And adjusting the fitness value of the updated particle by using the position and the speed of the updated particle swarm to obtain a global optimal value and an optimal solution.

In the process of finite element simulation of the three-dimensional ground stress field, firstly, the size of boundary load is preliminarily estimated from original field stress information obtained by a single-well rock mechanics model, then the size of the boundary load is gradually adjusted until an original field stress profile extracted by the three-dimensional ground stress model at a research well position is consistent with a single-well explained ground stress profile, and the extracted ground stress direction is ensured to be within certain precision with the ground stress trend obtained by each layer position imaging rose diagram formed by FMI imaging logging results of well drilling induced joints and well bore caving fracture trends at different positions.

The technical effect of the present invention is further described below in connection with a concrete simulation.

Simulation experiment

As shown in FIG. 3, the data diagram is applied to the natural fracture modeling results provided by the present invention. Fig. 4 is a graph of application data of rock mechanics numerical simulation considering natural fracture weakening provided by the present invention. FIG. 4(a) is a graph of the results of a three-dimensional rock force field model application without taking into account natural fracture weakening; fig. 4(b) is a graph of the results of a three-dimensional rock force field model application taking into account natural fracture weakening. Fig. 5 is a graph showing the application result of the three-dimensional ground stress field modeling provided by the present invention. FIG. 5(a) minimum horizontal principal stress distribution plot; (b) maximum horizontal principal stress profile. As shown in fig. 6, it is a comparison graph of the three-dimensional ground stress result and the one-dimensional ground stress magnitude result provided by the embodiment of the present invention.

Other embodiments of the disclosure will be apparent to those skilled in the art from consideration of the specification and practice of the disclosure disclosed herein. This application is intended to cover any variations, uses, or adaptations of the disclosure following, in general, the principles of the disclosure and including such departures from the present disclosure as come within known or customary practice within the art to which the disclosure pertains. It is intended that the specification and examples be considered as exemplary only, with a true scope and spirit of the disclosure being indicated by the following claims.

It will be understood that the present disclosure is not limited to the precise arrangements described above and shown in the drawings and that various modifications and changes may be made without departing from the scope thereof. The scope of the present disclosure should be limited only by the attached claims.

Claims (10)

1. A reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization is characterized by comprising the following steps of:

classifying the seismic facies according to the development characteristics of bedding cracks and natural cracks on the well;

under the constraint of the seismic facies, the crack lines are tracked by using the enhanced curvature attribute, and a deterministic model of a discrete natural crack network is established on the seismic scale;

based on the established crack certainty model, the grid conversion of the geological model and the ground stress finite element model is carried out by utilizing a related grid conversion program, and the optimal solution acquisition of the ground stress field is realized by utilizing a competitive particle swarm optimization algorithm.

2. The reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization according to claim 1, wherein the reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization specifically comprises the following steps:

step one, building a block structure model by using a well-seismic combination mode, comprising the following steps of: forming a whole-region construction model by using a two-dimensional seismic interpretation result as a transverse basis, using single-well layered data as a correction point and combining fault development characteristics;

determining longitudinal and transverse constraint conditions and selecting a variation function, and performing lithofacies modeling by using multi-point geostatistics and a sequential indication random simulation method; analyzing distribution ranges and frequency distribution forms of various attributes of different lithofacies on the basis of lithofacies modeling and structural modeling, performing key attribute modeling such as permeability, porosity and the like by an interpolation method, and performing constraint by single well interpretation data;

step three, utilizing a three-dimensional earthquake to describe and depict natural cracks of the block;

preprocessing a logging curve, correcting logging data of the expanding and well wall scouring positions, wherein the logging data comprise stratum gamma rays, density and longitudinal wave time, and rock parameters obtained through calculation according to acoustic logging data comprise elastic modulus and Poisson ratio key parameters; performing a uniaxial/triaxial compression test and a longitudinal and transverse sound wave speed test on the mechanical parameters of the rock core of the target layer, and analyzing and establishing a dynamic and static parameter conversion model;

step five, obtaining dynamic elastic body parameters through full-time window prestack three-dimensional seismic inversion from the earth surface to a target layer, calculating static data by utilizing a dynamic and static conversion model to obtain a three-dimensional rock mechanics data body, and introducing geological statistical information to specially predict thin layer elastic parameters for the thin layer mechanics parameters;

step six, introducing the fault serving as a discontinuous body into a three-dimensional rock mechanical parameter model, and analyzing the variation of the ground stress field at two sides of the fault; correcting the three-dimensional rock mechanical parameter model by combining single well logging, seismic interpretation results and indoor experimental results, and finally updating the three-dimensional rock mechanical field combined with well and seismic to obtain a three-dimensional rock mechanical parameter model aiming at geological structure influence;

step seven, utilizing single well acoustic-resistivity and density logging data to predict pore pressure, utilizing on-site drilling or logging actual measurement data to constrain, integrating drilling and well testing data, determining the interval between the maximum and minimum pore pressure values of a target layer, and constraining the pore pressure calculation result;

combining the characteristics of block geological stratification, reversely solving by using field collected small pressure data and well drilling data to obtain the structural stress coefficient of each layer, performing single-well logging ground stress calculation by using a combined spring model ground stress calculation formula, and correcting the logging calculation ground stress data by using the ground stress experimental results of different layers;

step nine, converting the geological corner point grid into a corresponding three-dimensional finite element grid model by using a grid interaction algorithm, searching node and unit information by using a global and local positioning method, and performing node assignment by using a spherical region searching method and a natural neighbor interpolation algorithm to realize the transmission of the geological model attribute result to a stress model of a finite element;

step ten, dispersing the block geometry;

step eleven, after determining a ground stress inversion work area and parameters of participation, establishing a three-dimensional ground stress field inversion target function, initializing particle swarm parameters, determining errors and iteration precision, and obtaining a global optimal solution by utilizing a competitive particle swarm algorithm;

and step twelve, carrying out finite element numerical simulation calculation, stopping calculation after iteration or calculation precision, and finally obtaining the three-dimensional ground stress distribution characteristics of the block.

3. The reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization according to claim 2, wherein the third step is to use a three-dimensional earthquake to describe and depict natural fractures of a block, and specifically comprises the following steps:

acquiring distribution characteristics of natural fractures by collecting multi-source information obtained by field outcrop, single-well core observation, logging cuttings and early geological research results in a research area, and constructing a bedding crack model by taking the established rock phase model and thickness reservoir model as constraints by combining with reservoir bedding weak plane explanation;

carrying out deterministic modeling on dispersion cracks by adopting a grouping mode to quantify the crack properties of the small-scale dispersion cracks;

performing topological parameters by adopting an ant colony tracking algorithm, establishing a related natural fracture network for large-scale discrete fracture analysis, forming a three-dimensional seismic data body by a plurality of two-dimensional seismic data surfaces, performing grid division on the two-dimensional seismic data surfaces and defining ant initial probability;

wherein, P_iIs a probability, C_iIs the coherence value at the ith point; under the condition of determining the starting probability, determining a main tracking direction by further dividing the region and calculating the gradient; the trace direction of each sub-block is expressed as

Wherein, theta_iTo track the direction G_x(i, j) and G_y(i, j) is the gradient of the data points.

4. The reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization according to claim 2, wherein the pre-stack elastic parameter inversion linear formula of the reflectivity of the step five based on direct inversion of Young's modulus, Poisson's ratio and density is as follows:

in the formula: r is a reflection coefficient; theta is an incident angle; v. of_PIs the velocity of the longitudinal wave; v. of_sIs the transverse wave velocity; k is

E is Young's modulus; sigma is the Poisson ratio; ρ is the density.

5. The reservoir three-dimensional stress field simulation method based on the deterministic multi-scale fracture modeling and the competitive particle swarm optimization algorithm according to claim 2, characterized in that in the sixth step, a mechanical parameter field after coupling a natural fracture is established aiming at the influence of fracture rigidity on rock mechanics;

Δμ_n: normal displacement, mm, k_n: normal stiffness, MPa/mm, k_s: shear stiffness, MPa/mm

Δτ_s: change in shear stress, MPa; Δ μ_s: tangential displacement, mm; delta tau_n: change in normal stress, MPa

The normal stress of the crack is in direct proportion to the rigidity of the crack, and the size of the crack is in inverse proportion to the rigidity of the crack; the rigidity parameter of the crack is drawn up according to the Bandis experimental data;

fracture shear stiffness is a function related to scale, fracture face normal stress, as follows:

in the formula: k_sShear stiffness of crack, σ_nThe normal stress of the crack surface, r the crack size, JRC the roughness, JCS the compressive strength,

is the rubbing angle.

6. The reservoir three-dimensional stress field simulation method based on the deterministic multi-scale fracture modeling and competitive particle swarm optimization algorithm according to claim 2, wherein in the seventh step, on-site drilling or logging actual measurement data are utilized for constraint, and if no blowout, well kick and description of obvious single gas occur in the drilling process, the drilling mud density is used as the upper limit of the formation pressure; and (3) balancing the stratum static pressure obtained in the well test with the stratum fluid in the shaft to obtain the stratum pressure, and if the balancing time is insufficient, taking the obtained static pressure as the lower limit of the stratum pressure.

7. The reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization according to claim 2, wherein in the eighth step, the combined spring model horizontal stress model calculation formula is as follows:

in the ninth step, a finite element load method is used for calculating a ground stress field, the three-dimensional research area is expanded in the horizontal direction, boundary grids are added, and the geometric size of each side of the expanded grid is 2-3 times that of the research area; meanwhile, a virtual grid is added, the three-dimensional grid extends upwards to the ground surface, and extends downwards to a proper depth below the reservoir stratum, so that the stress strain in the research area meets the continuous condition; meanwhile, introducing a new layer in the lowest layer of the stress model of the finite element, setting the bottom surface as a plane, and adding displacement constraint; taking the structure, the large fault and the fracture as stress model geometric boundaries of finite elements, and combining the gravity action to obtain the ground stress distribution of a horizontal stress boundary and a horizontal displacement boundary under the action of a load;

the step ten specifically comprises: the finite element grid adopts a random search method to carry out quality correction and treatment on the distorted and elongated grid so that the physical and mechanical properties of each rock are smooth and continuous between the grid and the grid, and finally, the amplitude of the parameter of the stratum material is obtained through well logging calculation data and well drilling actual measurement information.

8. The reservoir three-dimensional stress field simulation system based on the deterministic multi-scale fracture modeling and competitive particle swarm optimization based reservoir three-dimensional stress field simulation method according to any one of claims 1 to 7, is characterized by comprising the following steps:

the three-dimensional geological modeling module is used for obtaining related geological structures by interpreting seismic wave data of the whole bed system, correcting the seismic data by using single well data and geological stratification data, establishing an accurate block structure model by combining well-seismic information, modeling geological attributes by using a well-seismic combination modeling technology on the basis of three-dimensional structural modeling and phase modeling, and dividing fine layers of regions and platforms;

the natural fracture modeling module is used for constructing a bedding fracture model by taking the established rock phase model and the thickness reservoir model as constraints and combining with reservoir bedding weak plane explanation based on field outcrop, single-well core analysis, logging rock debris and early-stage geological analysis results of a research area; constructing a natural fracture model by combining the occurrence, density and development influence factors of the structural fractures in the research area, tracking fracture lines by using enhanced curvature attributes under the constraint of the seismic facies, realizing the establishment of a deterministic model of a discrete natural fracture network on the seismic scale, and correcting and debugging the natural fracture modeling result by combining near-wellbore geostatistics and well logging parameter fracture interpretation;

the logging curve correction and quality analysis module is used for performing single-well rock mechanical parameter and ground stress explanation by using a rock mechanical and ground stress model formula, and performing single-well geomechanical calculation result restriction by using relevant data such as a single-shaft/three-shaft compression test, a longitudinal and transverse sound wave speed test and the like of the rock mechanical parameter development;

the three-dimensional attribute model module is used for establishing a longitudinal and transverse wave and density three-dimensional attribute model of a research area by utilizing pre-stack seismic inversion, calculating the Young modulus and Poisson ratio of the rock, and obtaining a three-dimensional rock mechanical model on the basis of the distribution trend of seismic attribute bodies;

the three-dimensional rock mechanical field model correction module of the natural fracture is used for assigning measurement parameter values and model parameter initial values to different layers of the mechanical model according to the modeling results of the natural fracture and the fault and by combining the layer information of a single well; taking the fault as a discontinuous body, introducing the fault into a three-dimensional rock mechanical field model, and analyzing the change of the ground stress field at two sides of the fault; the method is also used for calculating the influence of the fracture rigidity on the rock mechanics parameter and establishing a three-dimensional rock mechanics parameter field after coupling the natural fracture;

the single-well stratum pore pressure model module is used for analyzing the deposition rule and the abnormal pressure generation mechanism of the research area, comprehensively utilizing logging data such as sound wave-resistivity and density to establish a single-well stratum pore pressure model, utilizing well drilling and completion and oil testing data to carry out data correction and quality control to obtain the relationship between the stratum pressure coefficient and the burial depth, lithology and porosity, and establishing the single-well stratum pore pressure model of the research area;

the grid conversion algorithm application module is used for searching node and unit information by adopting a global and local double positioning and tracking method, performing node assignment by adopting a spherical region searching method and a natural neighbor interpolation algorithm, converting geological corner grids into corresponding three-dimensional finite element grid models, and simultaneously utilizing a grid conversion reverse method and resampling the models to the corresponding corner grid geological models;

the rock stratum boundary unit adding module is used for adding overlying rock stratum, underlying rock stratum and lateral rock stratum boundary units for the research area, controlling the deformation of the reservoir stratum in a model space, and gradually reducing the displacement and stress concentration of the grid units along with the distance from the research area to reach the level of zero displacement and original site stress on the boundary;

the local ground stress field direction acquisition module is used for collecting drilling induced joints and well bore caving fracture trends of different layers of a work area, and forming an imaging rose diagram of each layer by combining FMI imaging logging results to obtain a local ground stress field direction;

the inversion grid boundary module is used for setting corresponding load and constraint conditions for a target reservoir, adding attributes for related media of the finite element model, performing elastoplasticity finite element calculation, and inverting the constraint, stress magnitude and direction conditions of the grid boundary;

and the competitive particle swarm algorithm solving module is used for acquiring a global optimal solution and a particle historical optimal solution by using the competitive particle swarm algorithm, comparing a calculation result with an actual measurement result of the ground stress size and direction of a known point, and stopping calculation after the accuracy or the iteration times are met.

9. An information data processing terminal, characterized in that the information data processing terminal comprises a memory and a processor, the memory stores a computer program, and the computer program when executed by the processor causes the processor to execute the reservoir three-dimensional stress field simulation method based on deterministic multi-scale fracture modeling and competitive particle swarm optimization according to any one of claims 1 to 7.

10. A computer readable storage medium storing a computer program which, when executed by a processor, causes the processor to execute the method for reservoir three-dimensional stress field simulation based on deterministic multi-scale fracture modeling and competitive particle swarm optimization as claimed in any one of claims 1 to 7.

Images