Elastic parameter brittleness evaluation numerical simulation method and device, electronic equipment and medium

Abstract

The invention discloses an elastic parameter brittleness evaluation numerical simulation method and device, electronic equipment and a medium. The method can comprise the steps that mesh generation is conducted on a shale model through a finite element method, and the constitutive relation of a solid matrix and a solid wave equation are obtained; obtaining a longitudinal wave speed calculation model, and calculating the longitudinal wave speed; obtaining a calculation model of the shear wave velocity, and calculating the shear wave velocity; obtaining a calculation model of the elastic parameters based on the longitudinal wave velocity and the transverse wave velocity; a brittleness index is calculated based on the elastic parameter. According to the invention, changes of elastic characteristics and brittleness characteristics of the shale with different bedding dip angles are researched in a numerical simulation mode, and influences of the elastic characteristics and brittleness indexes of the whole shale of the shale digital core with different bedding dip angle characteristics are researched through a digital rock physical experiment.

Description

Technical Field

[0001] The invention relates to the field of petroleum geophysical exploration and rock physics, and more specifically to a numerical simulation method, device, electronic equipment and medium for elastic parameter brittleness evaluation. Background Art

[0002] At present, some scholars have conducted research on hydraulic fracturing of shale reservoirs through a large number of physical experiments. Li Xiao et al. (2019) characterized the relationship between different shales' bedding inclinations and fracturing profiles through true triaxial fracturing experimental rock physics; through rock physics experiments, they studied the relationship between bedding's effects on the mechanical characteristics of rock, the failure mode and the bedding inclination. As the bedding inclination of shale increases, the uniaxial compressive strength and strain of shale show a trend of first decreasing and then increasing. Shale physics experiments are an important means of studying shale reservoir fracturing evaluation technology. However, when experimental rock physics is used to fractur e shale cores, the cores are fragile and non-repeatable, which limits the study of shale bedding dip and the overall elastic characteristics and brittleness index of shale. Therefore, the present invention uses the advantages of digital cores such as visibility and repeatability to carry out digital rock physics experimental research on shale bedding dip. Among them, the elastic parameter brittleness evaluation method is an important method for evaluating rock brittleness. It mainly uses Young's modulus and Poisson's ratio to calculate the brittleness index of rocks. Rocks with high Young's modulus and low Poisson's ratio usually have higher brittleness. The present invention applies it to digital rock physics experiments to obtain the brittleness characteristics of digital cores. Comprehensive research on shale bedding structure has not fully considered the influence of bedding dip on the overall elastic characteristics and brittleness index of shale, resulting in a lack of sufficient experimental support for the brittleness evaluation of shale reservoirs with different bedding dips.

[0003] At present, it is still necessary to develop a numerical simulation method for brittleness evaluation of shale elastic parameters with different bedding inclinations.

[0004] The information disclosed in the background technology section of the present invention is only intended to deepen the understanding of the general background technology of the present invention, and should not be regarded as acknowledging or suggesting in any form that the information constitutes the prior art already known to those skilled in the art. Summary of the invention

[0005] The present invention proposes a numerical simulation method, device, electronic equipment and medium for elastic parameter brittleness evaluation, which can study the changes in elastic characteristics and brittle characteristics of shales with different bedding inclinations through numerical simulation, and study the influence of the overall elastic characteristics and brittleness index of shale digital cores with different bedding inclination characteristics through digital rock physics experiments.

[0006] In a first aspect, an embodiment of the present disclosure provides a numerical simulation method for evaluating elastic parameter brittleness, comprising:

[0007] Finite element method is used to mesh the shale model, and the constitutive relation of the solid matrix and the solid wave equation are obtained;

[0008] Obtaining a calculation model for longitudinal wave velocity, and calculating the longitudinal wave velocity;

[0009] Obtaining a calculation model for shear wave velocity, and calculating the shear wave velocity;

[0010] Based on the longitudinal wave velocity and the shear wave velocity, a calculation model of elastic parameters is obtained;

[0011] Based on the elastic parameters, a brittleness index is calculated.

[0012] As a specific implementation of the embodiment of the present disclosure, the constitutive relationship of the solid matrix conforms to Hooke's law:

[0013] σ ij =λδ ij ε ii +2με ij

[0014] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

[0015] As a specific implementation of the embodiment of the present disclosure, the solid wave equation is:

[0016]

[0017] Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

[0018] As a specific implementation of the embodiment of the present disclosure, the calculation model of the longitudinal wave velocity is:

[0019]

[0020] Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

[0021] As a specific implementation of the embodiment of the present disclosure, the calculation model of the shear wave velocity is:

[0022]

[0023] Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

[0024] As a specific implementation of the embodiment of the present disclosure, the calculation model of the elastic parameter is:

[0025]

[0026]

[0027]

[0028] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0029] As a specific implementation of the embodiment of the present disclosure, the brittleness index is:

[0030]

[0031]

[0032]

[0033] Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

[0034] In a second aspect, the embodiments of the present disclosure further provide a numerical simulation device for evaluating elastic parameter brittleness, comprising:

[0035] The meshing module uses the finite element method to perform meshing on the shale model to obtain the constitutive relationship of the solid matrix and the solid wave equation;

[0036] A longitudinal wave velocity calculation module is used to obtain a calculation model of the longitudinal wave velocity and calculate the longitudinal wave velocity;

[0037] A shear wave velocity calculation module is used to obtain a shear wave velocity calculation model and calculate the shear wave velocity;

[0038] An elastic parameter calculation module, which obtains a calculation model of elastic parameters based on the longitudinal wave velocity and the transverse wave velocity;

[0039] The brittleness index calculation module calculates the brittleness index based on the elastic parameter.

[0040] As a specific implementation of the embodiment of the present disclosure, the constitutive relationship of the solid matrix conforms to Hooke's law:

[0041] σ ij =λδ ij ε ii +2με ij

[0042] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

[0043] As a specific implementation of the embodiment of the present disclosure, the solid wave equation is:

[0044]

[0045] Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

[0046] As a specific implementation of the embodiment of the present disclosure, the calculation model of the longitudinal wave velocity is:

[0047]

[0048] Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

[0049] As a specific implementation of the embodiment of the present disclosure, the calculation model of the shear wave velocity is:

[0050]

[0051] Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

[0052] As a specific implementation of the embodiment of the present disclosure, the calculation model of the elastic parameter is:

[0053]

[0054]

[0055]

[0056] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0057] As a specific implementation of the embodiment of the present disclosure, the brittleness index is:

[0058]

[0059]

[0060]

[0061] Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

[0062] In a third aspect, an embodiment of the present disclosure further provides an electronic device, the electronic device comprising:

[0063] A memory storing executable instructions;

[0064] A processor runs the executable instructions in the memory to implement the elastic parameter brittleness evaluation numerical simulation method.

[0065] In a fourth aspect, an embodiment of the present disclosure further provides a computer-readable storage medium, which stores a computer program, and when the computer program is executed by a processor, the numerical simulation method for evaluating elastic parameter brittleness is implemented.

[0066] Its beneficial effects are:

[0067] The present invention utilizes the advantages of digital cores such as visibility and repeatability to carry out digital rock physics experimental research on the relationship between shale bedding inclination and shale brittleness. It can well calculate the relationship between shale bedding structure and brittleness index, which is helpful to quantitatively carry out research on the relationship between shale bedding structure and shale reservoir elastic characteristics, and provides experimental support for shale reservoir brittleness evaluation.

[0068] The methods and apparatus of the present invention have other features and advantages that will be apparent from, or will be described in detail in, the accompanying drawings and subsequent detailed descriptions incorporated herein, which together serve to explain the specific principles of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0069] The above and other objects, features and advantages of the present invention will become more apparent through a more detailed description of exemplary embodiments of the present invention in conjunction with the accompanying drawings, wherein like reference numerals generally represent like components throughout the exemplary embodiments of the present invention.

[0070] Figure 1 A flow chart showing the steps of a numerical simulation method for elastic parameter brittleness evaluation according to an embodiment of the present invention.

[0071] Figure 2 A schematic diagram of the setting of longitudinal and transverse wave propagation boundary conditions according to an embodiment of the present invention is shown, wherein (a) represents the longitudinal wave propagation boundary condition, and (b) represents the transverse wave propagation boundary condition.

[0072] Figure 3 A schematic diagram showing the brittleness characteristics of different regions of a layered shale model according to an embodiment of the present invention.

[0073] Figure 4 A schematic diagram showing shale models with different bedding dips according to an embodiment of the present invention.

[0074] Figure 5 A schematic diagram showing the change of shale longitudinal wave velocity at different bedding angles according to an embodiment of the present invention is shown.

[0075] Figure 6 A schematic diagram showing the change of shale longitudinal wave velocity at different bedding angles according to an embodiment of the present invention is shown.

[0076] Figure 7 A schematic diagram showing the change of shale brittleness index at different bedding angles according to an embodiment of the present invention.

[0077] Figure 8A block diagram of a numerical simulation device for elastic parameter brittleness evaluation according to an embodiment of the present invention is shown.

[0078] Description of reference numerals:

[0079] 201. Grid generation module; 202. P-wave velocity calculation module; 203. S-wave velocity calculation module; 204. Elastic parameter calculation module; 205. Brittleness index calculation module. DETAILED DESCRIPTION

[0080] The preferred embodiments of the present invention will be described in more detail below. Although the preferred embodiments of the present invention are described below, it should be understood that the present invention can be implemented in various forms and should not be limited to the embodiments set forth herein.

[0081] To facilitate understanding of the solutions and effects of the embodiments of the present invention, six specific application examples are given below. Those skilled in the art should understand that the examples are only for facilitating understanding of the present invention, and any specific details thereof are not intended to limit the present invention in any way.

[0082] Example 1

[0083] Figure 1 A flow chart showing the steps of a numerical simulation method for elastic parameter brittleness evaluation according to an embodiment of the present invention.

[0084] like Figure 1 As shown, the numerical simulation method for elastic parameter brittleness evaluation includes: step 101, using the finite element method to perform grid division on the shale model to obtain the constitutive relationship of the solid matrix and the solid wave equation; step 102, obtaining the calculation model of the longitudinal wave velocity and calculating the longitudinal wave velocity; step 103, obtaining the calculation model of the shear wave velocity and calculating the shear wave velocity; step 104, based on the longitudinal wave velocity and the shear wave velocity, obtaining the calculation model of the elastic parameter; step 105, based on the elastic parameter, calculating the brittleness index.

[0085] In one example, the constitutive relation for a solid matrix is ​​Hooke's law:

[0086] σ ij =λδ ij ε ii +2με ij

[0087] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ijis the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

[0088] In one example, the wave equation for a solid is:

[0089]

[0090] Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

[0091] In one example, the calculation model of the longitudinal wave velocity is:

[0092]

[0093] Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

[0094] In one example, the shear wave velocity is calculated as:

[0095]

[0096] Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

[0097] In one example, the calculation model of the elasticity parameter is:

[0098]

[0099]

[0100]

[0101] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0102] In one example, the fragility index is:

[0103]

[0104]

[0105]

[0106] Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

[0107] Specifically, based on Hooke's law, solid wave equation and dynamic stress-strain method, the present invention obtains the changes in elastic parameters and brittleness index of shale models with different bedding angles, obtains the velocities of its longitudinal and transverse waves according to its stress-strain relationship, and further obtains the influence of bedding structure on the shale brittleness index.

[0108] The shale model is meshed using the finite element method, in which the constitutive relationship of the solid matrix in the shale conforms to Hooke's law:

[0109] σ ij =λδ ij ε ii +2με ij

[0110] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol, which is 1 when i=j and 0 when i≠j. The rock physics stress-strain experiment is simulated by the solid wave equation. The solid motion is given by the wave equation:

[0111]

[0112] Where ρ is the solid density, kg / m 3 ; u is the displacement field, m; S is the pressure field, Pa; F is the body force term.

[0113] By applying a sinusoidal displacement perpendicular to the top of the rock and parallel to the top of the rock, the stress-strain variation process of the layered shale is solved by numerical simulation, where the sinusoidal displacement conforms to:

[0114] u=A*sin(2πtf)

[0115] Where u is the displacement of the top of the rock, m; f is the frequency of the sinusoidal displacement, Hz; t is the time of the sinusoidal displacement, s; A is the amplitude of the sinusoidal displacement, which is 10 -9 m.

[0116] Figure 2 A schematic diagram showing the setting of boundary conditions for longitudinal and transverse wave propagation according to an embodiment of the present invention is shown.

[0117] When simulating the propagation of P- and S-waves in layered shales, the boundary condition setting is the key to affecting the propagation of P- and S-waves. The model boundary condition setting is the key to the effect of digital rock physics experiments in simulating the propagation of P- and S-waves. In order to meet the stability of the model, it is more in line with the conditions of digital rock physics experiments. Figure 2 As shown in the figure, when simulating the propagation of longitudinal waves, a sinusoidal vibration perpendicular to the top of the rock is applied to the top of the rock, where the X direction is the horizontal direction and the Y direction is the vertical direction. The average stress and average strain in the Y direction are calculated. The average strain in the Y direction is calculated by dividing the strain by the length of the rock in the Y direction. The effective P-wave modulus Re (M) is calculated from the average stress and strain at the top of the rock. The boundary solid displacement on both sides of the rock is 0 in the X direction, so that the strain can only occur in the longitudinal direction. The calculation model is as follows:

[0118]

[0119]

[0120] Where Re(M) represents the effective P-wave modulus, Pa; σ p represents the stress amplitude at the peak value during longitudinal wave propagation, m; ε p represents the strain amplitude at the peak value during longitudinal wave propagation, m; t 0 Indicates the time at the peak of the longitudinal wave propagation, s; V p is the longitudinal wave velocity, m / s; here we set the frequency f to 100 Hz.

[0121] When simulating shear wave propagation, a sinusoidal vibration parallel to the top of the rock is applied to the top of the rock, and the average stress and average strain in the X direction are calculated. The boundaries on both sides of the rock are periodic boundary conditions. When simulating longitudinal and transverse wave propagation, the displacement of the bottom boundary solid is 0. The calculation method of shear modulus Re (G) is similar to that of effective P-wave modulus, and both are calculated from the average stress-strain ratio of the top of the rock. The calculation model is as follows:

[0122]

[0123]

[0124] Where Re(G) represents the shear modulus, Pa; σ s represents the stress amplitude at the peak value during shear wave propagation, m; ε s represents the strain amplitude at the peak value during shear wave propagation, m; t 1Indicates the time at the peak of the shear wave propagation; V s is the shear wave velocity, m / s.

[0125] Through the above method, elastic parameters such as P- and S-wave velocities of shale models can be intuitively obtained through digital rock physics experiments, which provides a bridge for brittleness evaluation and joint analysis based on elastic parameters and mineral composition.

[0126] The dynamic Young's modulus, bulk modulus and Poisson's ratio can be calculated based on the longitudinal and transverse wave velocities obtained from the digital rock physics experiment, and the elastic parameters of the shale as a whole can be obtained. The calculation model is as follows:

[0127]

[0128]

[0129]

[0130] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0131] Using the elastic parameter brittleness evaluation method, the brittleness index B of shale is calculated after processing Young's modulus and Poisson's ratio. The calculation model is as follows:

[0132]

[0133]

[0134]

[0135] Where E is Young's modulus, MPSI; BRIT E Represents the effect of Young's modulus on the brittleness index; BRIT v Represents the effect of Poisson's ratio on the brittleness index. The range of Young's modulus is defined as 1MPSI to 8MPSI, and the range of Poisson's ratio is defined as 0.15 to 0.45. When E = 8MPSI and v = 0.15, the shale is considered to be "100%" brittle and has good fracturing properties; when E = 1MPSI and v = 0.45, the rock is considered to be "0%" brittle and has poor fracturing properties.

[0136] Example 2

[0137] The present invention also provides a numerical simulation device for evaluating elastic parameter brittleness, comprising:

[0138] The meshing module uses the finite element method to perform meshing on the shale model to obtain the constitutive relationship of the solid matrix and the solid wave equation;

[0139] A longitudinal wave velocity calculation module obtains a calculation model of the longitudinal wave velocity and calculates the longitudinal wave velocity;

[0140] A shear wave velocity calculation module is used to obtain a shear wave velocity calculation model and calculate the shear wave velocity;

[0141] The elastic parameter calculation module obtains the calculation model of elastic parameters based on the longitudinal wave velocity and the transverse wave velocity;

[0142] The brittleness index calculation module calculates the brittleness index based on the elastic parameters.

[0143] In one example, the constitutive relation for a solid matrix is ​​Hooke's law:

[0144] σ ij =λδ ij ε ii +2με ij

[0145] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

[0146] In one example, the wave equation for a solid is:

[0147]

[0148] Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

[0149] In one example, the calculation model of the longitudinal wave velocity is:

[0150]

[0151] Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

[0152] In one example, the shear wave velocity is calculated as:

[0153]

[0154] Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

[0155] In one example, the calculation model of the elasticity parameter is:

[0156]

[0157]

[0158]

[0159] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0160] In one example, the fragility index is:

[0161]

[0162]

[0163]

[0164] Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

[0165] Specifically, based on Hooke's law, solid wave equation and dynamic stress-strain method, the present invention obtains the changes in elastic parameters and brittleness index of shale models with different bedding angles, obtains the velocities of its longitudinal and transverse waves according to its stress-strain relationship, and further obtains the influence of bedding structure on the shale brittleness index.

[0166] The shale model is meshed using the finite element method, in which the constitutive relationship of the solid matrix in the shale conforms to Hooke's law:

[0167] σ ij =λδ ij ε ii +2με ij

[0168] Where λ and μ are the Lamé constants, σ ij and ε ijare the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol, which is 1 when i=j and 0 when i≠j. The rock physics stress-strain experiment is simulated by the solid wave equation. The solid motion is given by the wave equation:

[0169]

[0170] Where ρ is the solid density, kg / m 3 ; u is the displacement field, m; S is the pressure field, Pa; F is the body force term.

[0171] By applying a sinusoidal displacement perpendicular to the top of the rock and parallel to the top of the rock, the stress-strain variation process of the layered shale is solved by numerical simulation, where the sinusoidal displacement conforms to:

[0172] u=A*sin(2πtf)

[0173] Where u is the displacement of the top of the rock, m; f is the frequency of the sinusoidal displacement, Hz; t is the time of the sinusoidal displacement, s; A is the amplitude of the sinusoidal displacement, which is 10 -9 m.

[0174] When simulating the propagation of P- and S-waves in layered shales, the boundary condition setting is the key to affecting the propagation of P- and S-waves. The model boundary condition setting is the key to the effectiveness of digital rock physics experiments in simulating the propagation of P- and S-waves. In order to meet the stability of the model, it is more in line with the conditions of digital rock physics experiments.

[0175] Figure 2 A schematic diagram of the setting of longitudinal and transverse wave propagation boundary conditions according to an embodiment of the present invention is shown, wherein (a) represents the longitudinal wave propagation boundary condition, and (b) represents the transverse wave propagation boundary condition.

[0176] like Figure 2 As shown in (a), when simulating the propagation of longitudinal waves, a sinusoidal vibration perpendicular to the top of the rock is applied to the top of the rock, where the X direction is the horizontal direction and the Y direction is the vertical direction. The average stress and average strain in the Y direction are calculated. The average strain in the Y direction is calculated by dividing the strain by the length of the rock in the Y direction. The effective P-wave modulus Re (M) is calculated from the average stress and strain at the top of the rock. The boundary solid displacement on both sides of the rock is 0 in the X direction, so that the strain can only occur in the longitudinal direction. The calculation model is as follows:

[0177]

[0178]

[0179] Where Re(M) represents the effective P-wave modulus, Pa; σ p represents the stress amplitude at the peak value during longitudinal wave propagation, m; ε p represents the strain amplitude at the peak value during longitudinal wave propagation, m; t 0 Indicates the time at the peak of the longitudinal wave propagation, s; V p is the longitudinal wave velocity, m / s; here we set the frequency f to 100 Hz.

[0180] like Figure 2 As shown in (b), when simulating the propagation of shear waves, a sinusoidal vibration parallel to the top of the rock is applied to the top of the rock, and the average stress and average strain in the X direction are calculated. The boundaries on both sides of the rock are periodic boundary conditions. When simulating the propagation of longitudinal and shear waves, the displacement of the bottom boundary solid is 0. The calculation method of the shear modulus Re (G) is similar to that of the effective P-wave modulus, and both are calculated from the average stress-strain ratio of the top of the rock. The calculation model is as follows:

[0181]

[0182]

[0183] Where Re(G) represents the shear modulus, Pa; σ s represents the stress amplitude at the peak value during shear wave propagation, m; ε s represents the strain amplitude at the peak value during shear wave propagation, m; t 1 Indicates the time at the peak of the shear wave propagation; V s is the shear wave velocity, m / s.

[0184] Through the above method, elastic parameters such as P- and S-wave velocities of shale models can be intuitively obtained through digital rock physics experiments, which provides a bridge for brittleness evaluation and joint analysis based on elastic parameters and mineral composition.

[0185] The dynamic Young's modulus, bulk modulus and Poisson's ratio can be calculated based on the longitudinal and transverse wave velocities obtained from the digital rock physics experiment, and the elastic parameters of the shale as a whole can be obtained. The calculation model is as follows:

[0186]

[0187]

[0188]

[0189] Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0190] Using the elastic parameter brittleness evaluation method, the brittleness index B of shale is calculated after processing Young's modulus and Poisson's ratio. The calculation model is as follows:

[0191]

[0192]

[0193]

[0194] Where E is Young's modulus, MPSI; BRIT E Represents the effect of Young's modulus on the brittleness index; BRIT v Represents the effect of Poisson's ratio on the brittleness index. The range of Young's modulus is defined as 1MPSI to 8MPSI, and the range of Poisson's ratio is defined as 0.15 to 0.45. When E = 8MPSI and v = 0.15, the shale is considered to be "100%" brittle and has good fracturing properties; when E = 1MPSI and v = 0.45, the rock is considered to be "0%" brittle and has poor fracturing properties.

[0195] Example 3

[0196] Figure 3 A schematic diagram showing the brittleness characteristics of different regions of a layered shale model according to an embodiment of the present invention.

[0197] Construct a shale model with two layers, such as Figure 3 As shown, the model is a shale layered model with two bedding layers. The model is divided into three regions, among which the Young's modulus of region A and region C is 32.44 GPa and the Poisson's ratio is 0.23; the Young's modulus of region B is 35 GPa and the Poisson's ratio is 0.19. Compared with region B, region A and region C have lower Young's modulus and higher Poisson's ratio, showing a low brittle characteristic region, while region B has high Young's modulus and low Poisson's ratio, characterized by a high brittle characteristic region.

[0198] Figure 4 A schematic diagram showing shale models with different bedding dips according to an embodiment of the present invention.

[0199] The bedding inclination angles of shales are different, and the overall heterogeneity and anisotropy of shales are bound to be different. This paper studies the differences in the overall elastic characteristics of shales when the bedding is at different angles to the horizontal direction through digital rock physics experiments. Figure 4 As shown in the figure, a digital rock physics model of shale with different bedding dips was constructed. The contact relationship between the white and black areas of the model constitutes the bedding structure of the shale. The beddings are parallel to each other and the angles with the horizontal direction are 0°, 15°, 30°, 45°, 60°, 75°, and 90°. Figure 2The shale model in is two horizontal layers with an angle of 0 with the horizontal direction; Figure 4 The two layers of the shale model in (a) have an angle of 15° with the horizontal direction. The Young's modulus in the white area is 32.44 GPa and the Poisson's ratio is 0.23, and the Young's modulus in the red area is 35 GPa and the Poisson's ratio is 0.19. Similarly, Figure 4 In (b), (c), (d), (e), and (f), the angles between the two bedding layers and the horizontal direction of the shale model are 30°, 45°, 60°, 75°, and 90°, respectively. The Young's modulus in the white area is 32.44 GPa and the Poisson's ratio is 0.23, showing low brittleness characteristics. The Young's modulus in the black area is 35 GPa and the Poisson's ratio is 0.19, showing high brittleness characteristics.

[0200] Figure 5 A schematic diagram showing the change of shale longitudinal wave velocity at different bedding angles according to an embodiment of the present invention is shown.

[0201] Figure 6 A schematic diagram showing the change of shale longitudinal wave velocity at different bedding angles according to an embodiment of the present invention is shown.

[0202] Through the elastic parameter brittleness evaluation method, the P- and S-wave velocities and brittleness indexes of shale models with different bedding angles can be obtained, such as Figure 5 , Figure 6 shown. Figure 5 These are the P-wave and S-wave velocity diagrams of shale models with different bedding dip angles. It can be found that the P-wave velocity at dip angles of 60° and 75° is smaller than that of shale models with other bedding dip angles. The P-wave velocity of the shale model with a dip angle of 60° is 1934m / s, and the P-wave velocity of the shale with a dip angle of 75° is 1933m / s. When the dip angle is 90 degrees, the P-wave velocity of the model is roughly equal to the measured P-wave velocities of the models with dip angles of 0°, 15°, 30°, and 45°. From the relationship between the bedding dip and the shear wave velocity of the shale model, it can be seen that the shear wave velocity at dip angles of 60° and 75° is smaller than the shear wave velocity of the model at other dip angles. The shear wave velocity of the model at a dip angle of 60° is 1230m / s, and the shear wave velocity of the model at a dip angle of 75° is 1227m / s. The shear wave velocities measured by the model at dip angles of 0°, 15°, 30°, 45°, and 90° are roughly the same, and the values ​​are all around 1635m / s.

[0203] Figure 7 A schematic diagram showing the change of shale brittleness index at different bedding angles according to an embodiment of the present invention.

[0204] from Figure 7From the relationship between the dip angle of the middle bedding and the brittleness index of shale, it can be seen that when the dip angle is 60° and 75°, the brittleness index calculated by the model is less than the brittleness index calculated when the dip angle is 0°, 15°, 30°, and 45°. Since the difference in the P-wave and S-wave velocities of the shale model is small when the dip angle is 0°, 15°, 30°, and 45°, the brittleness index is roughly the same; when the dip angle is 60° and 75°, the P-wave and S-wave velocities calculated by the shale model are roughly the same, and the brittleness index calculated by these two bedding structures is also roughly the same. The results show that when shale contains the same material, the bedding dip angle is different, and the elastic characteristics and brittleness index of the whole shale will also be greatly different. When the dip angle of the shale bedding is 60° and 75° to the horizontal direction, the P-wave and S-wave velocities and brittleness index of the whole shale are compared with the shale model when the dip angle is 0°, 15°, 30°, and 45°. The calculated P-wave and S-wave velocities and brittleness index are smaller. Therefore, the different bedding directions of shale will have different effects on the brittleness index and related fracturing properties of shale. Therefore, the bedding inclination will also affect the fracturing properties of shale.

[0205] Example 4

[0206] Figure 8 A block diagram of a numerical simulation device for elastic parameter brittleness evaluation according to an embodiment of the present invention is shown.

[0207] like Figure 8 As shown, the elastic parameter brittleness evaluation numerical simulation device comprises:

[0208] A mesh generation module 201 is used to perform mesh generation on the shale model using the finite element method to obtain the constitutive relationship of the solid matrix and the solid wave equation;

[0209] A longitudinal wave velocity calculation module 202 obtains a calculation model of the longitudinal wave velocity and calculates the longitudinal wave velocity;

[0210] The shear wave velocity calculation module 203 obtains a shear wave velocity calculation model and calculates the shear wave velocity;

[0211] The elastic parameter calculation module 204 obtains a calculation model of the elastic parameters based on the longitudinal wave velocity and the shear wave velocity;

[0212] The brittleness index calculation module 205 calculates the brittleness index based on the elastic parameter.

[0213] As an alternative, the constitutive relation for the solid matrix is ​​Hooke's law:

[0214] σ ij =λδ ij ε ii +2με ij

[0215] Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

[0216] As an alternative, the solid wave equation is:

[0217]

[0218] Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

[0219] As an alternative, the calculation model of the longitudinal wave velocity is:

[0220]

[0221] Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

[0222] As an alternative, the calculation model of shear wave velocity is:

[0223]

[0224] Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

[0225] As an alternative, the calculation model of the elasticity parameter is:

[0226]

[0227]

[0228]

[0229] Among them, Kdyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

[0230] As an option, the brittleness index is:

[0231]

[0232]

[0233]

[0234] Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

[0235] Example 5

[0236] This embodiment provides an electronic device, which includes: a memory storing executable instructions; and a processor running the executable instructions in the memory to implement the above-mentioned elastic parameter brittleness evaluation numerical simulation method.

[0237] An electronic device according to an embodiment of the present disclosure includes a memory and a processor.

[0238] The memory is used to store non-temporary computer-readable instructions. Specifically, the memory may include one or more computer program products, which may include various forms of computer-readable storage media, such as volatile memory and / or non-volatile memory. The volatile memory may, for example, include random access memory (RAM) and / or cache memory (cache), etc. The non-volatile memory may, for example, include read-only memory (ROM), hard disk, flash memory, etc.

[0239] The processor may be a central processing unit (CPU) or other forms of processing units having data processing capabilities and / or instruction execution capabilities, and may control other components in the electronic device to perform desired functions. In one embodiment of the present disclosure, the processor is used to run the computer-readable instructions stored in the memory.

[0240] Those skilled in the art should be able to understand that in order to solve the technical problem of how to obtain a good user experience, the present embodiment may also include well-known structures such as a communication bus and an interface, and these well-known structures should also be included in the protection scope of the present disclosure.

[0241] For detailed description of this embodiment, reference may be made to the corresponding descriptions in the aforementioned embodiments, which will not be repeated here.

[0242] Example 6

[0243] This embodiment provides a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, the numerical simulation method for elastic parameter brittleness evaluation is implemented.

[0244] According to the computer-readable storage medium of the embodiment of the present disclosure, non-transitory computer-readable instructions are stored thereon. When the non-transitory computer-readable instructions are executed by a processor, all or part of the steps of the above-mentioned methods of each embodiment of the present disclosure are executed.

[0245] The above-mentioned computer-readable storage media include, but are not limited to: optical storage media (e.g., CD-ROM and DVD), magneto-optical storage media (e.g., MO), magnetic storage media (e.g., magnetic tape or mobile hard disk), media with built-in rewritable non-volatile memory (e.g., memory card) and media with built-in ROM (e.g., ROM box).

[0246] Those skilled in the art should understand that the purpose of the above description of the embodiments of the present invention is only to exemplarily illustrate the beneficial effects of the embodiments of the present invention, and is not intended to limit the embodiments of the present invention to any given examples.

[0247] The embodiments of the present invention have been described above, and the above description is exemplary, not exhaustive, and is not limited to the disclosed embodiments. Many modifications and changes will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments.

Claims

1. A numerical simulation method for evaluating brittleness of elastic parameters, It is characterized in that include: Finite element method is used to mesh the shale model, and the constitutive relation of the solid matrix and the solid wave equation are obtained; Obtaining a calculation model for longitudinal wave velocity, and calculating the longitudinal wave velocity; Obtaining a calculation model for shear wave velocity, and calculating the shear wave velocity; Based on the longitudinal wave velocity and the shear wave velocity, a calculation model of elastic parameters is obtained; Based on the elastic parameters, a brittleness index is calculated.

2. The numerical simulation method for elastic parameter brittleness evaluation according to claim 1, in, The constitutive relation of the solid matrix conforms to Hooke's law: s ij =λδ ij e ii +2with ij Where λ and μ are the Lamé constants, σ ij and ε ij are the elements in the stress tensor and strain tensor respectively, i and j are the positions of the elements in the tensor matrix, ε ii represents the elements in the strain tensor with the same number of rows and columns, δ ij is the Kronecker symbol. When i=j, its value is 1, and when i≠j, its value is 0.

3. The elastic parameter brittleness evaluation numerical simulation method according to claim 1, in, The solid wave equation is: Where ρ is the solid density, u is the displacement field, S is the pressure field, and F is the body force term.

4. The elastic parameter brittleness evaluation numerical simulation method according to claim 1, in, The calculation model of the longitudinal wave velocity is: Among them, V p is the longitudinal wave velocity, ρ is the solid density, Re(M) represents the effective longitudinal wave modulus, σ p represents the stress amplitude at the peak value during longitudinal wave propagation, ε p represents the strain amplitude at the peak value during longitudinal wave propagation, t 0 It represents the time at the peak value when the longitudinal wave propagates, and f is the frequency.

5. The elastic parameter brittleness evaluation numerical simulation method according to claim 1, in, The calculation model of the shear wave velocity is: Among them, V s is the shear wave velocity, ρ is the solid density, Re(G) represents the shear modulus, σ s represents the stress amplitude at the peak value during shear wave propagation, ε s represents the strain amplitude at the peak value during shear wave propagation, t 1 It represents the time at the peak of the shear wave propagation, and f is the frequency.

6. The elastic parameter brittleness evaluation numerical simulation method according to claim 1, in, The calculation model of the elastic parameters is: Among them, K dyn is the dynamic bulk modulus, E is the dynamic Young's modulus, and ν is the Poisson's ratio.

7. The elastic parameter brittleness evaluation numerical simulation method according to claim 1, in, The brittleness index is: Among them, B is the brittleness index, BRIT E Represents the influence of Young's modulus on the brittleness index, E is Young's modulus, BRIT v Represents the effect of Poisson's ratio on the brittleness index.

8. A numerical simulation device for evaluating elastic parameter brittleness, It is characterized in that include: The meshing module uses the finite element method to perform meshing on the shale model to obtain the constitutive relationship of the solid matrix and the solid wave equation; A longitudinal wave velocity calculation module is used to obtain a calculation model of the longitudinal wave velocity and calculate the longitudinal wave velocity; A shear wave velocity calculation module is used to obtain a shear wave velocity calculation model and calculate the shear wave velocity; An elastic parameter calculation module, which obtains a calculation model of elastic parameters based on the longitudinal wave velocity and the transverse wave velocity; The brittleness index calculation module calculates the brittleness index based on the elastic parameter.

9. An electronic device, It is characterized in that The electronic device comprises: A memory storing executable instructions; A processor, wherein the processor runs the executable instructions in the memory to implement the numerical simulation method for elastic parameter brittleness evaluation according to any one of claims 1 to 7.

10. A computer-readable storage medium, It is characterized in that The computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the numerical simulation method for evaluating elastic parameter brittleness according to any one of claims 1 to 7 is implemented.

Images