A refined prediction method for in-situ stress in shale formations under an observation coordinate system

Abstract

This invention discloses a refined prediction method for in-situ stress in shale formations under an observation coordinate system, comprising: obtaining basic physical property parameters through well logging data; establishing a shale rock physics model under the observation coordinate system; predicting and outputting the anisotropic stiffness coefficient and the verticalized P-wave and S-wave transit times under the constitutive coordinate system; calculating dynamic anisotropic rock mechanical parameters and static anisotropic rock mechanical parameters; calculating the formation pore pressure under the observation coordinate system; and combining the static anisotropic rock mechanical parameters and the predicted formation pore pressure to refine the evaluation of the maximum and minimum horizontal principal stresses. This invention establishes a refined prediction system for rock mechanical parameters, formation pore pressure, and in-situ stress in anisotropic shale formations under observation coordinate systems such as the build-up section and horizontal section, solving the problem that well logging data in anisotropic shale formations cannot be directly used for the evaluation of rock mechanical parameters, formation pore pressure, and in-situ stress under observation coordinate systems such as the build-up section and horizontal section.

Description

Technical Field

[0001] This invention relates to the field of shale formation stress technology, and in particular to a precise method for predicting shale formation stress under an observation coordinate system. Background Technology

[0002] Shale rocks are dense, with low porosity and permeability, and exhibit significant transverse isotropic characteristics (VTI). Shale gas reservoirs are typical artificial gas reservoirs, requiring hydraulic fracturing to achieve industrial production. To reduce costs, horizontal well drilling is currently the primary drilling technology for shale gas extraction. Creating hydraulic fractures through perforation and fracturing in the horizontal section is the main method for modifying shale formations. In-situ stress modeling of the horizontal section is crucial for precise fracturing design and optimization simulation. Single-well geomechanical modeling typically includes vertical stress calculation, formation pore pressure calculation, dynamic and static rock mechanics parameter calculation, and in-situ stress calculation. Among these, formation pore pressure prediction and dynamic and static rock mechanics parameter prediction are the key and challenging aspects of in-situ stress evaluation.

[0003] Predicting in-situ stress based on rock physics modeling is an important method. For example, Zhang Guangzhi et al. established an equivalent rock physics model for shale gas and used this model to calculate various elastic stiffness tensors, thereby predicting the minimum horizontal in-situ stress of shale formations. Liu Zhishui and Sun Zandong established a rock physics model of shale gas reservoirs considering kerogen particles and predicted the maximum and minimum horizontal principal stresses and fracture pressures of the reservoirs. Ma Tianshou et al. proposed to obtain shale physical property parameters, establish a transversely isotropic rock physics model, calculate stiffness coefficients, dynamic and static rock elastic parameters, and Biot coefficients, and use the Eaton method to predict formation pore pressure using P-wave transit time. Finally, they calculated the maximum and minimum horizontal principal stresses of transversely isotropic shale formations.

[0004] Because shale is characterized by low porosity, low permeability, and dense reservoirs, it requires hydraulic fracturing to achieve production capacity. Therefore, horizontal well fracturing is currently the main method for shale reservoir stimulation. The rock physics modeling methods for predicting in-situ stress mentioned above are all established in vertical wells. Due to the transverse isotropic nature of shale, the logging instrument and the symmetry axis of the shale formation no longer coincide during the logging process from the build-up section to the horizontal section. The sonic transit time measured at this time is the sonic transit time in the observation coordinate system, which differs from the sonic transit time measured in the constitutive coordinate system (logging instrument and shale formation symmetry axis) corresponding to vertical well logging. Therefore, current rock physics modeling methods for predicting in-situ stress are not applicable in observation coordinate systems such as horizontal wells. The sonic transit time measured in the observation coordinate system cannot be directly used for calculating dynamic anisotropic rock mechanics parameters and predicting formation pore pressure. Accurately evaluating the maximum and minimum horizontal principal stresses of shale formations in the observation coordinate system becomes even more difficult. Summary of the Invention

[0005] In view of this, the present invention provides a method for precise prediction of in-situ stress in shale formations under an observation coordinate system, so as to solve the above-mentioned technical problems.

[0006] This invention discloses a method for precise prediction of in-situ stress in shale formations under an observation coordinate system, which includes the following steps:

[0007] Step 1: Obtain basic physical property parameters through well logging data; wherein, the basic physical property parameters include acoustic wave, density, porosity, water saturation and shale mineral composition;

[0008] Step 2: Based on the basic physical property parameters, using the basic rock physics model, considering the observation coordinates of the wellbore trajectory, and combining the theory of anisotropic media, establish the shale rock physics model in the observation coordinate system;

[0009] Step 3: Based on the shale rock physics model in the observation coordinate system, predict and output the anisotropic stiffness coefficients and the longitudinal and transverse wave transit times in the constitutive coordinate system after verticalization.

[0010] Step 4: Based on the predicted anisotropic stiffness coefficient, calculate the dynamic anisotropic rock mechanical parameters, and combine the conversion relationship between the dynamic and static anisotropic rock mechanical parameters to calculate the static anisotropic rock mechanical parameters.

[0011] Step 5: Calculate the formation pore pressure in the observation coordinate system based on the P-wave time difference after verticalization.

[0012] Step 6: Combine static anisotropic rock mechanical parameters and predicted formation pore pressure to finely evaluate the maximum and minimum horizontal principal stresses.

[0013] Further, step 2 includes:

[0014] Step 21: Construct a physical model of brittle rock;

[0015] Step 22: Construct an equivalent elastic model of the dry clay-kerogen mixture;

[0016] Step 23: Construct a rock physics model for dry shale based on the brittle rock physics model and the equivalent elastic model of dry clay-kerogen mixture;

[0017] Step 24: Based on the gas saturation and water saturation, calculate the bulk modulus of the gas-water mixture in the pores using the Wood formula and convert it into a stiffness tensor; using the Brown-Korringa model, add the mixed fluid to the dry shale to obtain saturated fluid shale, and establish a transversely isotropic shale petrophysical model of saturated fluid shale.

[0018] Step 25: Based on the transversely isotropic shale rock physics model of the fluid, introduce the observation coordinate system and establish the shale rock physics model under the observation coordinate system.

[0019] Further, step 21 includes:

[0020] Step 211: Construct an equivalent elastic model for brittle mineral mixtures;

[0021] Step 212: Based on the equivalent elastic model of brittle mineral mixtures, establish a rock physical model of brittle shale.

[0022] Further, step 211 includes:

[0023] The formula for calculating the equivalent elastic modulus of brittle mineral mixtures using the Hashin-Shtrikman limit is as follows:

[0024]

[0025] In the formula, is the equivalent bulk modulus of a brittle mineral mixture, in GPa; The equivalent shear modulus for brittle mineral mixtures, GPa; , These are the maximum and minimum shear moduli in brittle minerals, expressed in GPa. , These represent the maximum and minimum bulk modulus values ​​for brittle minerals, expressed in GPa.

[0026] in, and Let be intermediate operators of the Hashingin-Shtrikman bound, and we have:

[0027]

[0028] In the formula, parentheses To calculate a weighted average of each brittle mineral based on its volume content;

[0029] Step 212 includes:

[0030] First, the equivalent elastic modulus of the dry shale matrix with an inorganic porosity of 50% was calculated using the isotropic SCA model:

[0031]

[0032] In the formula: ρ represents the bulk modulus of dry pores, in GPa; Here is the shear modulus of the dry pores, in GPa; is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of a brittle mineral mixture, in GPa; The equivalent bulk modulus of brittle rock with a porosity of 50% is given in GPa. The equivalent dry brittle rock shear modulus at a porosity of 50%, in GPa; , is the shape factor of the dry pores, which is dimensionless; , It is a shape factor of the matrix mineral mixture and is dimensionless.

[0033] Then, the absolute content of the brittle mineral mixture added to the background medium is: However, the relative content is The corresponding DEM model is:

[0034]

[0035] In the formula: The relative content of the remaining brittle mineral mixture to be added is dimensionless. The bulk modulus of brittle rock equivalent to that in the DEM model, in GPa; The equivalent dry rock shear modulus in the DEM model is given in GPa. is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of the matrix mineral mixture, in GPa; , The shape factor of the brittle mineral mixture as inclusions is dimensionless; the initial values ​​for iteration satisfy: , .

[0036] Further, step 22 includes:

[0037] Step 221: Simulate the equivalent elastic properties of clay-kerogen using the anisotropic SCA model and calculate the corresponding equivalent elastic stiffness matrix;

[0038] The calculation formula for the anisotropic SCA model is as follows:

[0039]

[0040] In the formula: Let GPa be the equivalent stiffness tensor for clay and kerogen. , , respectively, are the stiffness tensors of clay and kerogen, in GPa; It is a fourth-order unit tensor, dimensionless; , These are the Eshelby tensors for clay and kerogen, respectively, and are dimensionless. The relative content of clay relative to the clay-kerogen mixture;

[0041] Step 222: Establish a rock physics model of the dry clay-kerogen mixture;

[0042] Using a clay-kerogen mixture as the background medium, an anisotropic DEM model was established by adding empty organic pores to it. The corresponding anisotropic DEM model calculation formula is as follows:

[0043]

[0044] In the formula: The content of organic pores to be added to the background medium relative to the total amount of clay, kerogen, and organic pores is dimensionless. Let be the equivalent stiffness tensor of the dry clay-kerogen mixture, in GPa; Let GPa be the pore stiffness tensor. For the pore Eshelby stiffness tensor, which is dimensionless; The initial values ​​are fourth-order unit stiffness tensors, dimensionless; the initial values ​​satisfy the following during the iterative solution process: .

[0045] Further, step 23 includes:

[0046] Using Backus's average theory, a rock physics model of dry shale was established by mixing dry clay-kerogen mixtures with dry brittle rock, and the elastic stiffness coefficient of dry shale was calculated. The specific calculation formula is as follows:

[0047]

[0048] In the formula: The percentage of brittle rock to dry shale, a decimal, dimensionless; , , , , is the stiffness coefficient of dry shale, in GPa; , , , , is the stiffness coefficient of brittle rock, obtained from the establishment of a rock physical model of brittle shale, in GPa; , , , , is the stiffness coefficient of dry clay-kerogen, obtained from the establishment of the dry clay-kerogen model, in GPa.

[0049] Furthermore, in step 25:

[0050] The three wave velocities corresponding to the shale petrophysical model in the constructed observation coordinate system are expressed as follows:

[0051]

[0052] In the formula: c 11 , c 33 , c 44 , c 66 Here is the stiffness coefficient of shale rock in the constitutive coordinate system, in GPa; , , The P-wave velocity, SV-wave velocity, and SH-wave velocity of the shale rock physical model under the constructed observation coordinate system, in km / s; DEV The inclination angle is expressed in rad. ρ Density measured in well logging, g / cm³ 3 .

[0053] Furthermore, in step 3:

[0054] The formula for calculating the P-wave and S-wave time difference after verticalization is as follows:

[0055]

[0056] In the formula: c 33 , c 44 The shale rock stiffness coefficient in the constitutive coordinate system obtained for modeling, in GPa; DTCO_V and DTSM_V are the longitudinal and transverse wave transit times after verticalization, in us / ft; ρ Density measured in well logging, g / cm³ 3 .

[0057] Furthermore, in step 4:

[0058] The dynamic anisotropic rock mechanical parameters are calculated using the following formula:

[0059]

[0060] In the formula: The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and along the y-direction, and is dimensionless; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; , , , GPa represents the stiffness coefficient in the stiffness matrix.

[0061] By utilizing the conversion relationship between dynamic and static anisotropic Young's modulus, the static anisotropic rock mechanical parameters are calculated.

[0062] Further, step 5 includes:

[0063] Formation pore pressure is predicted using a formation pore pressure prediction method. The specific formation pore pressure prediction equation is as follows:

[0064]

[0065] In the formula, The predicted formation pore pressure is expressed in MPa. For vertical geostress, MPa; for pore pressure fitting coefficients of formations A and B, dimensionless; The verticalized P-wave time difference is expressed in µs / ft.

[0066] Step 6 includes:

[0067] The maximum and minimum horizontal principal stresses are precisely evaluated using the following calculation formula:

[0068]

[0069] In the formula, The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; The vertical ground stress is expressed in MPa. These are Biot coefficients, dimensionless; , These are the maximum and minimum strains in the horizontal direction, respectively, and are dimensionless. , The maximum and minimum horizontal principal stresses are predicted, in MPa.

[0070] Due to the adoption of the above technical solution, this invention has the following advantages: This invention establishes a refined prediction system for rock mechanical parameters, formation pore pressure, and in-situ stress in anisotropic shale formations under observation coordinate systems such as the build-up section and horizontal section. This solves the problem that well logging data in anisotropic shale formations cannot be directly used for the evaluation of rock mechanical parameters, formation pore pressure, and in-situ stress under observation coordinate systems such as the build-up section and horizontal section. The horizontal in-situ stress predicted by this invention under the observation coordinate system can provide important guidance for subsequent fracturing optimization design and in-situ fracturing operations. Attached Figure Description

[0071] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments recorded in the embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings.

[0072] Figure 1 This is a flowchart illustrating a method for fine prediction of in-situ stress in shale formations under an observation coordinate system, according to an embodiment of the present invention.

[0073] Figure 2 This is a schematic diagram of the physical property parameters of shale in well X1 according to an embodiment of the present invention;

[0074] Figure 3 This is a schematic diagram of a modeling method for anisotropic shale rock physics model in an observation coordinate system according to an embodiment of the present invention;

[0075] Figure 4 This is a schematic diagram of the stiffness coefficient and acoustic time difference verticalization inversion process according to an embodiment of the present invention;

[0076] Figure 5(a) is a schematic diagram of the horizontal to dynamic and static Young's modulus conversion relationship according to an embodiment of the present invention;

[0077] Figure 5(b) is a schematic diagram of the vertical dynamic and static Young's modulus conversion relationship according to an embodiment of the present invention;

[0078] Figure 6 This is a schematic diagram of the acoustic verticalization processing results and anisotropic rock mechanical parameter prediction results of an X1 well deviated well and horizontal well section according to an embodiment of the present invention.

[0079] Figure 7 This is a schematic diagram illustrating the prediction results of formation pore pressure and geostress in the X1 well deviated and horizontal well sections according to an embodiment of the present invention. Detailed Implementation

[0080] The present invention will be further described in conjunction with the accompanying drawings and embodiments. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of them. All other embodiments obtained by those skilled in the art should fall within the protection scope of the present invention.

[0081] See Figure 1 This invention provides an embodiment of a method for fine prediction of in-situ stress in shale formations under an observation coordinate system, which specifically includes:

[0082] S1. Using well logging data, obtain basic physical properties such as acoustic wave, density, porosity, water saturation, and shale mineral composition; see [link to article]. Figure 2 ;

[0083] S2. Based on the basic rock physics model, considering the observation coordinates of the wellbore trajectory, and combined with the theory of anisotropic media, a shale rock physics model under the observation coordinate system is constructed.

[0084] S21 Constructing a brittle rock physical model

[0085] The construction of a rock physics model for brittle shale involves two sub-steps.

[0086] S211 Constructs an equivalent elastic model for brittle mineral mixtures.

[0087] Considering that the elastic moduli of various brittle minerals do not differ significantly, this invention uses the Hashin-Shtrikman limit to simulate the upper and lower limits of the equivalent elastic modulus of brittle mineral mixtures, and takes the average of the upper and lower limits as the equivalent elastic modulus of the brittle mineral mixture. Compared with the Voigt-Reuss-Hill limit, the Hashin-Shtrikman limit has the narrowest upper and lower limits when simulating the equivalent elastic modulus of mineral mixtures; compared with the isotropic SCA model and the isotropic DEM model, the Hashin-Shtrikman limit does not need to consider the aspect ratio of mineral particles.

[0088] The formula for calculating the equivalent elastic modulus of brittle mineral mixtures using the Hashin-Shtrikman limit is shown below:

[0089]

[0090] In the formula: is the equivalent bulk modulus of a brittle mineral mixture, in GPa; The equivalent shear modulus for brittle mineral mixtures, GPa; , These are the maximum and minimum shear moduli in brittle minerals, expressed in GPa. , These are the maximum and minimum bulk modulus values ​​for brittle minerals, expressed in GPa.

[0091] and Let be intermediate operators of the Hashingin-Shtrikman bound, and we have:

[0092]

[0093] In the formula: parentheses To calculate a weighted average of each brittle mineral based on its volume content.

[0094] S212 Establishing a petrological model of brittle shale

[0095] Dry inorganic pores were added to a brittle mineral mixture using isotropic SCA and DEM models to obtain a petrophysical model of dry brittle shale. The combined use of the SCA and DEM models ensured the interconnectivity of the inorganic pores.

[0096] First, the equivalent elastic modulus of the dry shale matrix with an inorganic porosity of 50% was calculated using the isotropic SCA model.

[0097]

[0098] In the formula: ρ represents the bulk modulus of dry pores, in GPa; Here is the shear modulus of the dry pores, in GPa; is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of a brittle mineral mixture, in GPa; The equivalent bulk modulus of brittle rock with a porosity of 50% is given in GPa. The equivalent dry brittle rock shear modulus, in GPa, is given by a porosity of 50%. , is the shape factor of the dry pores, which is dimensionless; , It is a shape factor of the matrix mineral mixture and is dimensionless.

[0099] Then, brittle rock with a porosity of 50% was used as the background medium, and the porosity was adjusted to the actual state using a DEM model. Since the inorganic pores in brittle minerals are often less than 50%, in actual processing, it is necessary to reduce the porosity content and increase the content of the brittle mineral mixture; therefore, the additive is a matrix mineral mixture. It should be noted that the absolute content of the brittle mineral mixture added to the background medium is: However, the relative content should be Therefore, the corresponding DEM model is:

[0100]

[0101] In the formula: The relative content of the remaining brittle mineral mixture to be added is dimensionless. The bulk modulus of brittle rock equivalent to that in the DEM model, in GPa; The equivalent dry rock shear modulus in the DEM model is given in GPa. is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of the matrix mineral mixture, in GPa; , The shape factor of the brittle mineral mixture as inclusions is dimensionless; the initial values ​​for iteration satisfy: , .

[0102] S22 Constructs an equivalent elastic model of dry clay-kerogen mixture.

[0103] It includes two sub-steps:

[0104] S221 calculates the corresponding equivalent elastic stiffness matrix.

[0105] The equivalent elastic properties of clay-kerogen were simulated using an anisotropic SCA model, and the corresponding equivalent elastic stiffness matrix was calculated.

[0106] The calculation formula for the anisotropic SCA model is as follows:

[0107]

[0108] In the formula: Let GPa be the equivalent stiffness tensor for clay and kerogen. , , respectively, are the stiffness tensors of clay and kerogen, in GPa; It is a fourth-order unit tensor, dimensionless; , These are the Eshelby tensors for clay and kerogen, respectively, and are dimensionless. This represents the relative content of clay relative to the clay-kerogen mixture.

[0109] Unlike the method proposed by Gui Junchuan et al. (2020) to simulate the elastic properties of clay-kerogen by combining anisotropic SCA model and anisotropic DEM model, this paper argues that using only the anisotropic SCA model is sufficient to simulate the equivalent elastic characteristics of clay-kerogen mixtures.

[0110] S222 Establishing a rock physics model of dry clay-kerogen mixture

[0111] Using a clay-kerogen mixture as the background medium, an anisotropic DEM model was established by adding empty organic pores to it. The corresponding anisotropic DEM model calculation formula is as follows:

[0112]

[0113] In the formula: The content of organic pores to be added to the background medium relative to the total amount of clay, kerogen, and organic pores is dimensionless. Let be the equivalent stiffness tensor of the dry clay-kerogen mixture, in GPa; Let GPa be the pore stiffness tensor. For the pore Eshelby stiffness tensor, which is dimensionless; The initial values ​​are fourth-order unit stiffness tensors, dimensionless; the initial values ​​satisfy the following during the iterative solution process: .

[0114] Construction of the petrological model of S23 dry shale

[0115] Using Backus's average theory, a rock physics model of dry shale was established by mixing dry clay-kerogen mixtures with dry brittle rocks, and the elastic stiffness coefficient of dry shale was calculated. The specific calculation formula is shown below:

[0116]

[0117] In the formula: The percentage of brittle rock to dry shale, a decimal, dimensionless; , , , , is the stiffness coefficient of dry shale, in GPa; , , , , is the stiffness coefficient of brittle rock, obtained from the establishment of a rock physical model of brittle shale, in GPa; , , , , is the stiffness coefficient of dry clay-kerogen, obtained from the establishment of the dry clay-kerogen model, in GPa.

[0118] Using the Backus model to mix brittle rock and dry clay-kerogen avoids the problem of needing to select one or more minerals as the background medium.

[0119] Establishment of a shale rock physical model in the S24 constitutive coordinate system

[0120] S241 Calculate the bulk modulus of a gas-water mixture.

[0121] Based on the gas saturation and water saturation, the bulk modulus of the gas-water mixture in the pores is calculated using the Wood formula (Wood, 1955) and then converted into a stiffness tensor.

[0122] Wood's formula (Wood, 1955) is shown below:

[0123]

[0124] In the formula: Let be the bulk modulus of the gas, in GPa; Let be the bulk modulus of water, in GPa; Let be the equivalent bulk modulus of the mixed fluid, in GPa; The gas saturation is dimensionless. Let GPa be the equivalent bulk modulus of the mixed fluid.

[0125] S242 Fluid Replacement

[0126] Considering the transversely isotropic nature of dry shale, the Brown-Korringa model is used to add mixed fluids to dry shale to obtain saturated fluid shale, and a petrophysical model of transversely isotropic shale with saturated fluids is established.

[0127] The calculation formula for the Brown-Korringa model is as follows:

[0128]

[0129] In the formula: For the compliance tensor of saturated fluid shale, GPa -1 ; For the flexibility tensor of dry shale, GPa -1 ; For the compliance tensor of the shale matrix, GPa -1 ; , The compressibility coefficient, distributed as fluid and rock matrix, in GPa -1 ; Porosity is dimensionless.

[0130] Establishment of a shale rock physical model in the S25 observation coordinate system

[0131] Introducing an observation coordinate system and constructing a shale petrological model within that coordinate system, see [link to relevant documentation]. Figure 3In VTI shale formations, wellbore inclination angle and formation dip angle are important factors leading to differences in the elastic properties of shale between the observation coordinate system and the constitutive coordinate system. Assuming the shale formation is gently sloping (dip angle of 0°), the wellbore inclination angle DEV is the angle between the wellbore axis and the shale bedding symmetry axis in the observation coordinate system. In this case, the three wave velocities corresponding to the constructed shale petrophysical model in the observation coordinate system can be expressed as:

[0132]

[0133] In the formula: c 11 , c 33 , c 44 , c 66 Here is the stiffness coefficient of shale rock in the constitutive coordinate system, in GPa; , , The P-wave velocity, SV-wave velocity, and SH-wave velocity of the shale rock physical model under the constructed observation coordinate system, in km / s; DEV The inclination angle is expressed in rad. ρ Density measured in well logging, g / cm³ 3 .

[0134] The P-wave, fast S-wave, and slow S-wave measured in the well logging coordinate system in the build-up and horizontal sections correspond to the P-wave, SH-wave, and SV-wave mentioned above, respectively.

[0135] The shale mineral components and fluid elastic moduli required in the modeling process are shown in Table 1.

[0136] Table 1 Elastic parameters of major minerals and fluids in shale (Mavko et al., 2009; Carcione et al., 2011)

[0137]

[0138] S3. Based on the shale rock physics model in the observation coordinate system, and using the sonic transit time of well logging as a constraint, simulated annealing is employed to predict and output the stiffness coefficients and vertically adjusted P-wave and S-wave transit times in the constitutive coordinate system. (See also...) Figure 4 .

[0139] The formulas for calculating the P-wave and S-wave time difference after verticalization are as follows:

[0140]

[0141] In the formula: c 33 , c44 The shale rock stiffness coefficient in the constitutive coordinate system obtained for modeling, in GPa; DTCO_V and DTSM_V are the longitudinal and transverse wave transit times after verticalization, in us / ft; ρ Density measured in well logging, g / cm³ 3 .

[0142] S4. Calculate the dynamic anisotropic rock mechanics parameters using the predicted anisotropic stiffness coefficients, and calculate the static anisotropic rock mechanics parameters by combining the transformation relationship between dynamic and static anisotropic rock mechanics parameters.

[0143]

[0144] In the formula: The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and along the y-direction, and is dimensionless; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; , , , GPa represents the stiffness coefficient in the stiffness matrix.

[0145] Considering that the dynamic and static Poisson's ratios differ in magnitude, no dynamic-to-static conversion is performed on the Poisson's ratio. The static anisotropic rock mechanical parameters are calculated using the conversion relationship between the dynamic and static anisotropic Young's moduli. See Figures 5(a) and 5(b).

[0146] S5. Calculate the formation pore pressure in the observation coordinate system using the longitudinal wave time difference after verticalization.

[0147] Formation pore pressure is predicted using a formation pore pressure prediction method. The specific formation pore pressure prediction equation is as follows:

[0148]

[0149] In the formula, The predicted formation pore pressure is expressed in MPa. For vertical geostress, MPa; for pore pressure fitting coefficients of formations A and B, dimensionless; The verticalized P-wave time difference is expressed in µs / ft.

[0150] S6. Based on anisotropic static rock mechanics parameters and predicted formation pore pressure, a detailed evaluation of the maximum and minimum horizontal principal stresses is conducted. (See also...) Figure 6 .

[0151] The method for precise evaluation of in-situ stress is as follows:

[0152]

[0153] In the formula, The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; The vertical ground stress is expressed in MPa. These are Biot coefficients, dimensionless; , These are the maximum and minimum strains in the horizontal direction, respectively, and are dimensionless. , The maximum and minimum horizontal principal stresses are predicted, in MPa.

[0154] The present invention also provides a specific embodiment:

[0155] Using well logging data, shale physical properties parameters are obtained, including: shale mineral composition, water saturation, density, porosity, and P-wave and S-wave transit times in the observation coordinate system, such as... Figure 2 As shown.

[0156] Using a basic rock physics model, and fully considering the anisotropic origin of shale, a shale rock physics model in a constitutive coordinate system is constructed. Based on this model, and considering the influence of the wellbore trajectory observation coordinate system on the logging sonic transit time, an anisotropic shale rock physics model in an observation coordinate system is constructed by introducing the theory of anisotropic shale media. Figure 3 As shown.

[0157] by Figure 2 The shale physical properties shown are the basic input data, and the data is used to... Figure 3 A shale rock physics model in the constructed observation coordinate system predicts the acoustic transit time in the observation coordinate system, and based on simulated annealing, utilizes... Figure 4 The stiffness coefficient and acoustic transit time verticalization inversion process shown involves repeatedly comparing and iterating the predicted acoustic transit time in the observation coordinate system with the acoustic transit time measured in the well logging until the error between the two meets the accuracy requirements. The process then outputs five stiffness coefficients corresponding to the constitutive coordinate system. c 11 , c 13 , c 33 , c 44 , c 66And the longitudinal wave time difference DTCO_V after verticalization treatment, such as Figure 6 As shown in lanes 6 and 4.

[0158] Using the stiffness coefficients in the obtained constitutive coordinate system, the dynamic horizontal and vertical Young's moduli and Poisson's ratio are calculated. Then, based on the dynamic-static Young's modulus conversion relationship shown in Figure 5(a) and Figure 5(a) respectively, the static horizontal and vertical Young's moduli are calculated. Considering the poor correlation and mutual magnitude between the dynamic and static Poisson's ratios, no dynamic-static conversion is performed on the Poisson's ratio. The predicted static anisotropic rock mechanical parameters are as follows: Figure 6 As shown in the 7th passage.

[0159] Vertical in-situ stress was calculated using density logging data, and based on the verticalized P-wave transit time DTCO_V, formation pore pressure in the observation coordinate system was predicted using empirical formulas for predicting formation pore pressure based on vertical well sections within the work area. Figure 7 As shown in the 7th question.

[0160] Based on the calculated horizontal and vertical Young's modulus and Poisson's ratio, combined with the vertical in-situ stress and predicted formation pore pressure, the maximum and minimum horizontal principal stresses in the observation coordinate systems of the build-up section and horizontal sections are precisely evaluated using the maximum and minimum horizontal principal stress calculation method established based on the vertical well section within the work area. Figure 7 As shown in the 8th question.

[0161] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit it. Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the specific implementation of the present invention. Any modifications or equivalent substitutions that do not depart from the spirit and scope of the present invention should be covered within the scope of protection of the claims of the present invention.

Claims

1. A method for precise prediction of in-situ stress in shale formations under an observation coordinate system, characterized in that, Includes the following steps: Step 1: Obtain basic physical property parameters through well logging data; wherein, the basic physical property parameters include acoustic wave, density, porosity, water saturation and shale mineral composition; Step 2: Based on the basic physical property parameters, using the basic rock physics model, considering the observation coordinates of the wellbore trajectory, and combining the theory of anisotropic media, establish the shale rock physics model in the observation coordinate system; Step 3: Based on the shale rock physics model in the observation coordinate system, predict and output the anisotropic stiffness coefficients and the longitudinal and transverse wave transit times in the constitutive coordinate system after verticalization. Step 4: Based on the predicted anisotropic stiffness coefficient, calculate the dynamic anisotropic rock mechanical parameters, and combine the conversion relationship between the dynamic and static anisotropic rock mechanical parameters to calculate the static anisotropic rock mechanical parameters. Step 5: Calculate the formation pore pressure in the observation coordinate system based on the P-wave time difference after verticalization. Step 6: Combine static anisotropic rock mechanical parameters and predicted formation pore pressure to finely evaluate the maximum and minimum horizontal principal stresses; Step 2 includes: Step 21: Construct a physical model of brittle rock; Step 22: Construct an equivalent elastic model of the dry clay-kerogen mixture; Step 23: Construct a rock physics model for dry shale based on the brittle rock physics model and the equivalent elastic model of dry clay-kerogen mixture; Step 24: Based on the gas saturation and water saturation, calculate the bulk modulus of the gas-water mixture in the pores using the Wood formula and convert it into a stiffness tensor; using the Brown-Korringa model, add the mixed fluid to the dry shale to obtain saturated fluid shale, and establish a transversely isotropic shale petrophysical model of saturated fluid shale. Step 25: Based on the transversely isotropic shale rock physics model of the fluid, introduce the observation coordinate system and establish the shale rock physics model under the observation coordinate system; In step 25: The three wave velocities corresponding to the shale petrophysical model in the constructed observation coordinate system are expressed as follows: In the formula: c 11 , c 33 , c 44 , c 66 Here is the stiffness coefficient of shale rock in the constitutive coordinate system, in GPa; , , The P-wave velocity, SV-wave velocity, and SH-wave velocity of the shale rock physical model under the constructed observation coordinate system are given in km / s. DEV The inclination angle is expressed in rad. ρ Density measured in well logging, g / cm³ 3 ; In step 3: The formula for calculating the P-wave and S-wave time difference after verticalization is as follows: In the formula: c 33 , c 44 The shale rock stiffness coefficient in the constitutive coordinate system obtained for modeling, in GPa; DTCO_V and DTSM_V are the longitudinal and transverse wave transit times after verticalization, in us / ft; ρ Density measured in well logging, g / cm³ 3 ; Step 5 includes: Formation pore pressure is predicted using a formation pore pressure prediction method. The specific formation pore pressure prediction equation is as follows: In the formula, The predicted formation pore pressure is expressed in MPa. For vertical geostress, MPa; for pore pressure fitting coefficients of formations A and B, dimensionless; The verticalized P-wave time difference is expressed in µs / ft.

2. The method according to claim 1, characterized in that, Step 21 includes: Step 211: Construct an equivalent elastic model for brittle mineral mixtures; Step 212: Based on the equivalent elastic model of brittle mineral mixtures, establish a rock physical model of brittle shale.

3. The method according to claim 2, characterized in that, Step 211 includes: The formula for calculating the equivalent elastic modulus of brittle mineral mixtures using the Hashin-Shtrikman limit is as follows: In the formula, is the equivalent bulk modulus of a brittle mineral mixture, in GPa; The equivalent shear modulus for brittle mineral mixtures, GPa; , These are the maximum and minimum shear moduli in brittle minerals, expressed in GPa. , These represent the maximum and minimum bulk modulus values ​​for brittle minerals, expressed in GPa. in, and Let be intermediate operators of the Hashingin-Shtrikman bound, and we have: In the formula, parentheses To calculate a weighted average of each brittle mineral based on its volume content; Step 212 includes: First, the equivalent elastic modulus of the dry shale matrix with an inorganic porosity of 50% was calculated using the isotropic SCA model: In the formula: ρ represents the bulk modulus of dry pores, in GPa; Here is the shear modulus of the dry pores, in GPa; is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of a brittle mineral mixture, in GPa; The equivalent bulk modulus of brittle rock with a porosity of 50% is given in GPa. The equivalent dry brittle rock shear modulus at a porosity of 50%, in GPa; , is the shape factor of the dry pores, which is dimensionless; , It is a shape factor of the matrix mineral mixture and is dimensionless. Then, the absolute content of the brittle mineral mixture added to the background medium is: However, the relative content is The corresponding DEM model is: In the formula: The relative content of the remaining brittle mineral mixture to be added is dimensionless. The bulk modulus of brittle rock equivalent to that in the DEM model, in GPa; The equivalent dry rock shear modulus in the DEM model, in GPa; is the equivalent bulk modulus of a brittle mineral mixture, in GPa; is the equivalent shear modulus of the matrix mineral mixture, in GPa; , The shape factor of the brittle mineral mixture as inclusions is dimensionless; the initial values ​​for iteration satisfy: , .

4. The method according to claim 1, characterized in that, Step 22 includes: Step 221: Simulate the equivalent elastic properties of clay-kerogen using the anisotropic SCA model and calculate the corresponding equivalent elastic stiffness matrix; The calculation formula for the anisotropic SCA model is as follows: In the formula: Let GPa be the equivalent stiffness tensor for clay and kerogen. , , respectively, are the stiffness tensors of clay and kerogen, in GPa; It is a fourth-order unit tensor, dimensionless; , These are the Eshelby tensors for clay and kerogen, respectively, and are dimensionless. The relative content of clay relative to the clay-kerogen mixture; Step 222: Establish a rock physics model of the dry clay-kerogen mixture; Using a clay-kerogen mixture as the background medium, an anisotropic DEM model was established by adding empty organic pores to it. The corresponding anisotropic DEM model calculation formula is as follows: In the formula: The content of organic pores to be added to the background medium relative to the total amount of clay, kerogen, and organic pores is dimensionless. Let be the equivalent stiffness tensor of the dry clay-kerogen mixture, in GPa; Let GPa be the pore stiffness tensor. For the pore Eshelby stiffness tensor, which is dimensionless; The initial values ​​are fourth-order unit stiffness tensors, dimensionless; the initial values ​​satisfy the following during the iterative solution process: .

5. The method according to claim 1, characterized in that, Step 23 includes: Using Backus's average theory, a rock physics model of dry shale was established by mixing dry clay-kerogen mixtures with dry brittle rock, and the elastic stiffness coefficient of dry shale was calculated. The specific calculation formula is as follows: In the formula: The percentage of brittle rock to dry shale, a decimal, dimensionless; , , , , ρ is the stiffness coefficient of dry shale, in GPa; , , , , is the stiffness coefficient of brittle rock, obtained from the establishment of a rock physical model of brittle shale, in GPa; , , , , is the stiffness coefficient of dry clay-kerogen, obtained from the establishment of the dry clay-kerogen model, in GPa.

6. The method according to claim 1, characterized in that, In step 4: The dynamic anisotropic rock mechanical parameters are calculated using the following formula: In the formula: The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and along the y-direction, and is dimensionless; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; , , , GPa represents the stiffness coefficient in the stiffness matrix. By utilizing the conversion relationship between dynamic and static anisotropic Young's modulus, the static anisotropic rock mechanical parameters are calculated.

7. The method according to claim 1, characterized in that, Step 6 includes: The maximum and minimum horizontal principal stresses are precisely evaluated using the following calculation formula: In the formula, The Young's modulus parallel to the bedding plane is GPa; The Young's modulus perpendicular to the bedding plane, GPa; Poisson's ratio is parallel to the bedding plane and has no dimensions; Poisson's ratio perpendicular to the bedding plane, dimensionless; The vertical ground stress is expressed in MPa. These are Biot coefficients, dimensionless; , These are the maximum and minimum strains in the horizontal direction, respectively, and are dimensionless. , The maximum and minimum horizontal principal stresses are predicted, in MPa.

Images