A method for in-situ stress inversion using MCMC and stress polygon constraints

Abstract

The present invention belongs to the technical field of deep formation drilling and completion, and discloses a method for geostress inversion using MCMC and stress polygon constraints, including collecting prior distribution information of geomechanical parameters, well wall damage information, and observed uniaxial compressive strength of rocks and their standard deviations; sampling the target geomechanical parameters based on the proposed distribution under the condition of satisfying the quadratic constraint model of geostress based on the MCMC method; calculating the acceptance rate of the sampling results based on the prior distribution information and the observed samples using a delayed rejection algorithm; updating the proposed distribution characteristics using an adaptive algorithm; repeating the above sampling process until the sampling number is met, and discarding the combustion period to achieve geostress inversion. The present invention takes into account the influence of uncertainty of reservoir geomechanical parameters, and more objectively analyzes the probability distribution of horizontal geostress in the reservoir through the combination of prior information and observed samples, providing theoretical and technical guidance for the prediction of horizontal geostress and wellbore stress analysis.

Description

Technical Field

[0001] The invention relates to a geostress inversion method using MCMC and stress polygon constraints, and belongs to the technical field of oil and gas drilling and completion. Background Art

[0002] Geostress is mainly the force per unit area generated inside the medium caused by the movement of the earth's crust and the complex underground environment. It is usually characterized by three principal stresses: overlying rock pressure, horizontal maximum geostress and horizontal minimum geostress. For (ultra) deep formations, the effective determination of geostress will directly affect the evaluation of wellbore stability, formation sand analysis and the implementation of hydraulic fracturing operations.

[0003] Among them, the overburden rock pressure can be easily predicted numerically by the rock density integration method; while the prediction methods of horizontal geostress mainly include indoor core testing technology, single well testing and empirical data analysis and other measurement technologies and analysis methods. Specific field measurement methods include paleomagnetic core guidance technology, acoustic emission Kaiser effect, small hydraulic fracturing and well collapse analysis. However, in engineering, it is more common to use statistics and fitting based on theoretical calculation formulas to predict horizontal geostress. In addition, the inherent component differences of the formation rock itself and the uncertainty influence in the sampling and testing process will always lead to the emergence of geomechanical parameter uncertainty, which further reduces the prediction accuracy of the existing theoretical horizontal geostress calculation. In this regard, it has become a common understanding to analyze underground engineering problems by incorporating uncertainties such as inherent fuzziness and experimental errors of geomechanical parameters. However, most of the existing theoretical prediction models of geostress ignore the characteristics of parameter uncertainty, and do not organically combine historical prior uncertainty with observation samples.

[0004] At present, there are many theoretical calculation methods for predicting geostress. For example, the invention patent with publication number CN116990862A of CNOOC (China) Co., Ltd. Shenzhen Branch discloses a method for predicting geostress under the constraint of rock mechanical layers. This method abandons the disadvantages of large differences in stress field distribution caused by single lithology and interlayer lithology mutations, and focuses on the interlayer lithology and its interlayer bonding relationship. By dividing the rock mechanical layers and calculating the equivalent mechanical parameters, the geostress prediction under complex lithology background is realized.

[0005] Although the above-mentioned patents and other existing technical solutions have achieved certain technical effects, they ignore the uncertainty characteristics of parameters, resulting in a decrease in the accuracy of ground stress prediction in complex strata, thereby affecting the implementation of underground projects. Summary of the invention

[0006] In order to overcome the defects existing in the prior art, the present invention aims to provide a geostress inversion method using MCMC and stress polygon constraints. The present invention can achieve an effective combination of geostress prior information and observation data, so that the predicted geostress can be given in the form of probabilistic uncertainty taking into account uncertainty, thereby solving the problem of geostress prediction accuracy containing geomechanical parameter uncertainty.

[0007] The technical solution provided by the present invention to solve the above technical problems is: a method for inversion of geostress using MCMC and stress polygon constraints, comprising the following steps:

[0008] According to the historical information of the block, logging data and indoor experimental results, the prior distribution information of geomechanical parameters, well wall damage information, observed uniaxial compressive strength of rock and standard deviation are collected;

[0009] A quadratic constraint model of ground stress for wellbore damage information is established based on stress polygon and wellbore stability evaluation model;

[0010] Based on the MCMC method, the target geomechanical parameters are sampled according to the proposed distribution under the condition of satisfying the quadratic constraint model of ground stress;

[0011] The delayed rejection algorithm is used to calculate the acceptance rate of the sampling results according to the prior distribution information and the observed samples, and the sampling chain is updated;

[0012] Adopting adaptive algorithm to update the proposed distribution features;

[0013] The sampling process is repeated until the sampling number is met, and the combustion period is discarded to achieve the inversion of ground stress.

[0014] A further technical solution is that the prior distribution information of geomechanical parameters, well wall damage information and observed uniaxial compressive strength of rock and its standard deviation collected based on block historical information, logging data and indoor experimental results include:

[0015] Collect prior distribution information of geomechanical parameters based on historical information of the block;

[0016] Obtain well wall damage information based on imaging logging technology and engineering data;

[0017] The geomechanical parameters of the wellbore damage area, the uniaxial compressive strength of rock and its standard deviation are obtained using logging data and indoor experimental results.

[0018] A further technical solution is that the calculation formula of the standard deviation is:

[0019]

[0020] Where: For the iUniaxial compressive strength of samples, MPa; is the average value of uniaxial compressive strength, MPa; is the sample size; is the standard deviation.

[0021] A further technical solution is that the quadratic constraint model of geostress includes a stress polygon constraint geostress model and a wellbore stability evaluation model.

[0022] A further technical solution is that the stress polygon constrained geostress model is:

[0023]

[0024] Where: is the friction coefficient with a value between 0.6 and 1.0, dimensionless; is the maximum horizontal ground stress, MPa; is the minimum horizontal ground stress, MPa; is the overlying rock pressure, MPa; is the pore pressure, MPa.

[0025] A further technical solution is that the wellbore stability evaluation model includes:

[0026] Rock stress field around the well:

[0027]

[0028] Where: , , , , , is the effective stress of the rock around the well in cylindrical coordinate system, MPa; , , , , , is the in-situ stress field variable in the rectangular coordinate system at the wellbore, MPa; is the drilling fluid column pressure, MPa; is the pore pressure, MPa; v is Poisson’s ratio, dimensionless; is the well circumference failure angle, °;

[0029] Rock shear failure criterion:

[0030]

[0031] Where: is the minimum effective principal stress of the rock around the well, MPa; is the maximum effective principal stress of the rock around the well, MPa; is the uniaxial compressive strength, MPa; is the internal friction angle, °;

[0032] Tensile rupture criteria:

[0033]

[0034] Where: S t is the tensile strength of rock, MPa, is the minimum effective principal stress of the rock around the well, MPa.

[0035] A further technical solution is that the sampling of the target geomechanical parameters based on the MCMC method according to the proposed distribution under the condition of satisfying the quadratic constraint model of ground stress includes:

[0036] The proposed distribution is used to perform Monte Carlo sampling on the target geomechanical parameters under the condition of satisfying the quadratic constraint model of geostress, and a stable Markov is constructed by observing the samples. The sampling results approach the posterior distribution of the target parameters; the proposed distribution refers to t The transfer kernel of the Markov chain at time is in the form of satisfying the mean , the covariance is Normal distribution .

[0037] A further technical solution is that the acceptance rate calculation formula is:

[0038]

[0039]

[0040] Where: Likelihood function is the prior distribution, For the i Observed samples of rock strength, MPa; The rock strength predicted by the wellbore stability evaluation model, MPa; It is the sampling result of the last moment; For this moment in the recommended distribution The sampling results below are is the first-stage acceptance rate; u is the random number sampled in [0,1] in the first stage at this moment; is the second stage acceptance rate; u 2 is the random number sampled in [0,1] in the second stage at this moment.

[0041] A further technical solution is that the use of an adaptive algorithm to update the suggested distribution features includes:

[0042] Current running times t Less than the set number of times t When 0, the proposed distribution maintains the initial covariance C 0 unchanged; when the number of runs t > t 0, the proposed distribution updates the covariance C t for:

[0043]

[0044] Where: s d is dependent on the dimension of the parameter space d The constant of C 0 is the initial covariance; t is the number of runs; t 0 is the set number of times; C t is the covariance; For the i Sampling results of rock strength.

[0045] A further technical solution is that the inversion of the ground stress by discarding the combustion period is specifically: removing the unstable samples in the initial stage of the MCMC sampling process to obtain the ground stress inversion result.

[0046] The present invention also provides an electronic device, comprising:

[0047] one or more processors;

[0048] a storage device having one or more programs stored thereon;

[0049] When the one or more programs are executed by the one or more processors, the one or more processors are enabled to implement the in-situ stress inversion method using MCMC and stress polygon constraints as described above in the present invention.

[0050] The present invention also provides a storage medium having a computer program stored thereon, wherein when the computer program is executed by a processor, the inversion method of geostress using MCMC and stress polygon constraints as described above is implemented in the present invention.

[0051] The present invention has the following beneficial effects: the inversion method of geostress using MCMC and stress polygon constraints provided by the present invention establishes a geostress quadratic constraint model about wellbore damage information based on stress polygons and a wellbore stability evaluation model, wherein the MCMC method is used to achieve effective fusion of geomechanical parameter prior distribution information and observation data, and the delayed rejection (DR) algorithm and the adaptive (AM) algorithm are used to improve the sampling quality and efficiency; thereby, the present invention can automatically sample the posterior distribution of geostress under the condition of parameter uncertainty, providing theoretical support for predicting the numerical value of deep formation stress. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] Figure 1 It is a sampling chain of horizontal ground stress obtained by MCMC inversion in an embodiment of the present invention;

[0053] Figure 2 It is the quantitative result of the uncertainty of horizontal ground stress obtained by inverting the wellbore damage information in the embodiment of the present invention. DETAILED DESCRIPTION

[0054] The technical solution of the present invention will be described clearly and completely below in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of them. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0055] A method for in-situ stress inversion using MCMC and stress polygon constraints of the present invention comprises the following steps:

[0056] Step 1: Collect prior distribution information of geomechanical parameters based on historical information of the block;

[0057] The prior distribution information of geomechanical parameters mainly includes: using the historical analysis data of the block to obtain the mechanical parameters of the formation rock (including elastic modulus and Poisson's ratio), rock strength parameters (including compressive strength, internal friction angle and tensile strength) and the probability distribution range of the ground stress and pore pressure values;

[0058] Step 2: Obtain well wall damage information based on imaging logging technology and engineering data;

[0059] The wellbore damage information refers to the use of imaging logging technology to observe and determine the type of damage at the wellbore, namely, collapse / rupture, and obtain information such as the damage angle and collapse width at that location based on the imaging logging technology;

[0060] Step 3: Use logging data and indoor experimental results to obtain the geomechanical parameters of the wellbore damage area and the uniaxial compressive strength of the rock and its standard deviation;

[0061]

[0062] Where: For the i Uniaxial compressive strength of samples, MPa; is the average value of uniaxial compressive strength, MPa; is the sample size; is the standard deviation;

[0063] Step 4: Establish a geostress quadratic constraint model for wellbore damage information based on the stress polygon and wellbore stability evaluation model;

[0064] The quadratic constrained geostress model includes the stress polygon constrained geostress model and the wellbore stability evaluation model that uses the rock stress field around the well to substitute the rock failure criteria (shear failure and tensile rupture):

[0065] Among them, the stress polygon constrained ground stress model is:

[0066]

[0067] Where: is the friction coefficient with a value between 0.6 and 1.0, dimensionless; is the maximum horizontal ground stress, MPa; is the minimum horizontal ground stress, MPa; is the overlying rock pressure, MPa; is the pore pressure, MPa;

[0068] The rock stress field around the well:

[0069]

[0070] Where: , , , , , is the effective stress of the rock around the well in cylindrical coordinate system, MPa; , , , , , is the in-situ stress field variable in the rectangular coordinate system at the wellbore, MPa; is the drilling fluid column pressure, MPa; is the pore pressure, MPa; v is Poisson’s ratio, dimensionless; is the well circumference failure angle, °;

[0071] The failure angle around the well is taken as , ,

[0072]

[0073] Where: is the initial collapse angle, °, is the avalanche width, °; is the initial failure angle of shear fracture, °;

[0074] for , The value of is assumed that the uniaxial compressive strength and tensile strength of rock are unknown variables. After uniformly taking values ​​at certain intervals within the wellbore angle [0°, 360°], the shear failure criterion and tensile failure criterion are respectively combined in the rock stress field around the well to obtain the rock uniaxial compressive strength matrix and tensile strength matrix required for the wellbore wall to maintain stability. The wellbore angle at which the maximum uniaxial compressive strength and tensile strength are required corresponds to , and the maximum uniaxial compressive strength and tensile strength are the minimum constraint requirements for maintaining the stability of the wellbore under the current ground stress state.

[0075] , , , , , is the in-situ stress field variable in the rectangular coordinate system at the wellbore, which is obtained by the following relationship:

[0076]

[0077] Where: is the maximum horizontal ground stress, MPa; is the minimum horizontal ground stress, MPa; is the overlying rock pressure, MPa; w is the angle between the maximum geostress and the wellbore azimuth, °; is the well inclination angle;

[0078] The rock shear failure criterion:

[0079]

[0080]

[0081] Where: is the minimum effective principal stress of the rock around the well, MPa; is the maximum effective principal stress of the rock around the well, MPa; is the uniaxial compressive strength, MPa; is the internal friction angle, °; C 0 is cohesion, MPa;

[0082] The tensile rupture criterion:

[0083]

[0084] Where: S t is the tensile strength of rock, MPa, is the minimum effective principal stress of the rock around the well, MPa.

[0085] The calculation formula for the maximum effective principal stress and minimum effective principal stress of the rock around the well is:

[0086]

[0087] Step 5: Based on the MCMC method, sample the target geomechanical parameters according to the current recommended distribution under the condition of satisfying the quadratic constraint model of geostress;

[0088] The MCMC method is to combine the Markov chain with Monte Carlo simulation, construct a Markov process by observing samples, so that the stationary distribution in the simulated Markov chain approaches the parameter posterior distribution. The theoretical form of the Markov chain is:

[0089]

[0090] Where: X N is a random sequence, where ; is the state parameter set;

[0091] Step 6: Calculate the sampling acceptance rate based on the prior information and the observed sample based on the MCMC method combined with the delayed rejection (DR) algorithm. If the sample is judged to be accepted, update the sampling chain ;

[0092] The formula for calculating the acceptance rate is:

[0093]

[0094]

[0095] Where: Likelihood function is the prior distribution, For the i Observed samples of rock strength, MPa; The rock strength predicted by the wellbore stability evaluation model, MPa; is the sampling result of the previous moment, and at this moment is equal to the mean of the recommended distribution described in step 5 ; For this moment in the recommended distribution The sampling results below are is the first-stage acceptance rate; u is the random number sampled in [0,1] in the first stage at this moment; is the second stage acceptance rate; u 2 is the random number sampled in [0,1] in the second stage at this moment;

[0096] Step 7: Adopt the adaptive (AM) algorithm to update the current proposed distribution features;

[0097] Specifically: Current running times t Less than the set number of times t When 0, the proposed distribution maintains the initial covariance C 0 unchanged; when the number of runs t > t 0, the proposed distribution updates the covariance C t for:

[0098]

[0099]

[0100] Where: s d is dependent on the dimension of the parameter space d The constant of C 0 is the initial covariance; t is the number of runs; t 0 is the set number of times; C t is the covariance; For the i Sampling results of rock strength;

[0101] Step 8: Repeat the sampling process of steps 5 to 7 until the sampling number is met. According to the MCMC sampling chain, remove the unstable samples in the initial stage of the MCMC sampling process to obtain the inversion results of the geostress.

[0102] Example

[0103] In this embodiment, the parameters involved in the in-situ stress inversion method using MCMC and stress polygon constraints (the relevant parameters are derived from the X well 4860.50 m where the wellbore collapse occurred but no rupture occurred) are shown in Table 1:

[0104] Table 1

[0105]

[0106] In this embodiment, the maximum horizontal principal stress is set to and the horizontal minimum principal stress is the target geomechanical parameter, and its distribution characteristic is a non-informative prior distribution.

[0107] According to the above parameters, we can get Figure 1 The MCMC sampling results for the burner are shown, Figure 2 The uncertainty quantification results of horizontal ground stress are obtained by inverting in the ground stress quadratic constraint model excluding the burning period.

[0108] The embodiment of the present invention further provides a geostress inversion system using MCMC and stress polygon constraints, which is used to implement the geostress inversion method using MCMC and stress polygon constraints as described above, including:

[0109] Geostress quadratic constraint model calculation module: used to calculate the current stress field around the well and evaluate the wellbore stability; and analyze whether the sampling results of geostress meet the stress polygon constraint conditions and the actual wellbore damage information under the current wellbore damage information and rock strength conditions;

[0110] MCMC-based sampling module: used to calculate the acceptance rate based on prior information and observed samples using the delayed rejection (DR) algorithm, and to update the current proposed distribution characteristics using the adaptive (AM) algorithm, and finally obtain the inversion results of the ground stress after removing the burning period.

[0111] The embodiments of the present invention also provide corresponding electronic devices and computer-readable storage media for implementing the solutions provided by the embodiments of the present invention.

[0112] The electronic device includes a storage device and a processor, the storage device is used to store instructions or codes, and the processor is used to execute the instructions or codes so that the electronic device executes the geostress inversion method using MCMC and stress polygon constraints described in any embodiment of the present application.

[0113] In practical applications, the computer-readable storage medium may be any combination of one or more computer-readable media, which may be a computer-readable signal medium or a computer-readable storage medium.

[0114] The above description is not intended to limit the present invention in any form. Although the present invention has been disclosed through the above embodiments, it is not intended to limit the present invention. Any technician familiar with the profession can make some changes or modifications to equivalent embodiments of equivalent changes using the technical contents disclosed above without departing from the scope of the technical solution of the present invention. However, any simple modification, equivalent changes and modifications made to the above embodiments based on the technical essence of the present invention without departing from the content of the technical solution of the present invention are within the scope of the technical solution of the present invention.

Claims

1. A geostress inversion method using MCMC and stress polygon constraints, characterized in that: The following steps are involved: According to the historical information of the block, logging data and indoor experimental results, the prior distribution information of geomechanical parameters, well wall damage information, observed uniaxial compressive strength of rock and standard deviation are collected; A quadratic constraint model of ground stress for wellbore damage information is established based on stress polygon and wellbore stability evaluation model; Based on the MCMC method, the target geomechanical parameters are sampled according to the proposed distribution under the condition of satisfying the quadratic constraint model of ground stress; The proposed distribution is used to perform Monte Carlo sampling on the target geomechanical parameters under the condition of satisfying the quadratic constraint model of geostress, and a stable Markov is constructed by observing the samples. The sampling results approach the posterior distribution of the target parameters; the proposed distribution refers to t The transfer kernel of the Markov chain at time is in the form of satisfying the mean , the covariance is Normal distribution ; The delayed rejection algorithm is used to calculate the acceptance rate of the sampling results according to the prior distribution information and the observed samples, and the sampling chain is updated; The acceptance rate calculation formula is: Where: Likelihood function is the prior distribution, For the i Observed samples of rock strength, MPa; The rock strength predicted by the wellbore stability evaluation model, MPa; It is the sampling result of the last moment; For this moment in the recommended distribution The sampling results below are is the first-stage acceptance rate; u is the random number sampled in [0,1] in the first stage at this moment; is the second stage acceptance rate; u 2 is the random number sampled in [0,1] in the second stage at this moment; Adopting adaptive algorithm to update the proposed distribution features; The sampling process is repeated until the sampling number is met, and the combustion period is discarded to achieve the inversion of ground stress.

2. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 1, characterized in that: The prior distribution information of geomechanical parameters, well wall damage information and observed uniaxial compressive strength of rock and its standard deviation collected based on historical information of the block, logging data and indoor experimental results include: Collect prior distribution information of geomechanical parameters based on historical information of the block; Obtain well wall damage information based on imaging logging technology and engineering data; The geomechanical parameters of the wellbore damage area, the uniaxial compressive strength of rock and its standard deviation are obtained using logging data and indoor experimental results.

3. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 2, characterized in that: The calculation formula of the standard deviation is: Where: For the i Uniaxial compressive strength of samples, MPa; is the average value of uniaxial compressive strength, MPa; is the sample size; is the standard deviation.

4. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 1, characterized in that: The geostress quadratic constraint model includes a stress polygon constraint geostress model and a wellbore stability evaluation model.

5. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 4, characterized in that: The stress polygon constrained geostress model is: Where: is the friction coefficient with a value between 0.6 and 1.0, dimensionless; is the maximum horizontal ground stress, MPa; is the minimum horizontal ground stress, MPa; is the overlying rock pressure, MPa; is the pore pressure, MPa.

6. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 4, characterized in that: The wellbore stability evaluation model includes: Rock stress field around the well: Where: , , , , , is the effective stress of the rock around the well in cylindrical coordinate system, MPa; , , , , , is the in-situ stress field variable in the rectangular coordinate system at the wellbore, MPa; is the drilling fluid column pressure, MPa; is the pore pressure, MPa; v is Poisson's ratio, dimensionless; is the well circumference failure angle, °; Rock shear failure criterion: Where: is the minimum effective principal stress of the rock around the well, MPa; is the maximum effective principal stress of the rock around the well, MPa; is the uniaxial compressive strength, MPa; is the internal friction angle, °; Tensile rupture criteria: Where: S t is the tensile strength of rock, MPa, is the minimum effective principal stress of the rock around the well, MPa.

7. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 1, characterized in that: The use of an adaptive algorithm to update the suggested distribution features includes: Current running times t Less than the set number of times t When 0, the proposed distribution maintains the initial covariance C 0 unchanged; when the number of runs t > t 0, the proposed distribution updates the covariance C t for: Where: s d is dependent on the dimension of the parameter space d The constant of C 0 is the initial covariance; t is the number of runs; t 0 is the set number of times; C t is the covariance; For the i Sampling results of rock strength.

8. The in-situ stress inversion method using MCMC and stress polygon constraints according to claim 1, characterized in that: The inversion of the ground stress by discarding the burning period is specifically as follows: the unstable samples in the initial stage of the MCMC sampling process are removed to obtain the ground stress inversion result.

Images