Seismic inversion method and device, storage medium and electronic equipment

Abstract

This invention provides a seismic inversion method, apparatus, storage medium, and electronic device. The method includes: acquiring initial inversion data, and acquiring first azimuth seismic data and second azimuth seismic data; constructing at least one first azimuth forward modeling operator and at least one second azimuth forward modeling operator; constructing a first target inversion function based on at least one first azimuth forward modeling operator and the first azimuth seismic data; and constructing a second target inversion function based on at least one second azimuth forward modeling operator and the second azimuth seismic data; determining the target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function, and the initial inversion data; and calculating the target inversion result of the anisotropy parameters based on the target inversion result of the parameters to be inverted. Embodiments of this invention can improve the accuracy of seismic inversion.

Description

Technical Field

[0001] This invention relates to the field of seismic exploration of oil and gas, and particularly to a seismic inversion method, apparatus, storage medium and electronic equipment. Background Technology

[0002] Seismic waves in subsurface media typically exhibit azimuth anisotropy of varying types and intensities, with reflection and transmission responses also varying with azimuth. Horizontal Transverse Isotropy (HTI) media, as a medium model inducing azimuth anisotropy, is one of the most effective anisotropy models and plays a crucial role in describing seismic wave propagation. Anisotropy is significant for describing fracture orientation and intensity, as well as the type of filling fluid. However, fracture-induced anisotropy contributes relatively little to seismic reflection, leading to significant ambiguity in predictions directly from azimuth seismic data. To address this, related techniques using inverted stiffness matrices to calculate anisotropy parameters suffer from ill-posedness, resulting in low accuracy in seismic inversion. Therefore, a satisfactory solution for improving the accuracy of seismic inversion remains elusive. Summary of the Invention

[0003] In view of this, embodiments of the present invention provide a seismic inversion method, apparatus, storage medium, and electronic device to solve the problem of low accuracy in seismic inversion using related technologies. That is, embodiments of the present invention can achieve seismic inversion through a first objective inversion function and a second objective inversion function to predict the anisotropy parameters of the HTI medium, thereby obtaining a highly accurate objective inversion result of the anisotropy parameters, which can effectively improve the accuracy of seismic inversion.

[0004] According to one aspect of the present invention, a seismic inversion method is provided, the method comprising:

[0005] Acquire initial inversion data, as well as first azimuth seismic data and second azimuth seismic data;

[0006] Construct at least one first azimuth forward operator and at least one second azimuth forward operator;

[0007] Based on the at least one first azimuth forward modeling operator and the first azimuth seismic data, a first objective inversion function is constructed; and based on the at least one second azimuth forward modeling operator and the second azimuth seismic data, a second objective inversion function is constructed.

[0008] Based on the first target inversion function, the second target inversion function, and the initial inversion data, the target inversion result of the parameter to be inverted is determined;

[0009] Based on the target inversion results of the parameters to be inverted, the target inversion results of the anisotropy parameters are calculated.

[0010] According to another aspect of the present invention, a seismic inversion apparatus is provided, the apparatus comprising:

[0011] The acquisition unit is used to acquire initial inversion data, as well as first azimuth seismic data and second azimuth seismic data;

[0012] The processing unit is used to construct at least one first azimuth angle forward operator and at least one second azimuth angle forward operator;

[0013] The processing unit is further configured to construct a first objective inversion function based on the at least one first azimuth forward modeling operator and the first azimuth seismic data; and to construct a second objective inversion function based on the at least one second azimuth forward modeling operator and the second azimuth seismic data;

[0014] The processing unit is further configured to determine the target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function, and the initial inversion data;

[0015] The processing unit is further configured to calculate the target inversion result of the anisotropy parameters based on the target inversion result of the parameters to be inverted.

[0016] According to another aspect of the present invention, an electronic device is provided, the electronic device including a processor and a memory storing a program, wherein the program includes instructions that, when executed by the processor, cause the processor to perform the methods mentioned above.

[0017] According to another aspect of the present invention, a non-transitory computer-readable storage medium storing computer instructions for causing a computer to perform the methods mentioned above is provided.

[0018] According to another aspect of the present invention, a computer program product is provided, comprising a computer program, wherein the computer program, when executed by a processor, is used to cause the computer to perform the methods mentioned above.

[0019] This invention, after acquiring initial inversion data and first and second azimuth seismic data, constructs at least one first azimuth forward modeling operator and at least one second azimuth forward modeling operator. Then, based on at least one first azimuth forward modeling operator and the first azimuth seismic data, a first target inversion function is constructed; and based on at least one second azimuth forward modeling operator and the second azimuth seismic data, a second target inversion function is constructed. Further, based on the first and second target inversion functions and the initial inversion data, the target inversion results for the parameters to be inverted are determined; and based on the target inversion results for the parameters to be inverted, the target inversion results for anisotropy parameters are calculated. Therefore, this invention can achieve seismic inversion through the first and second target inversion functions to predict the anisotropy parameters of HTI media, thereby obtaining highly accurate target inversion results for anisotropy parameters and effectively improving the accuracy of seismic inversion. Attached Figure Description

[0020] Further details, features, and advantages of the invention are disclosed in the following description of exemplary embodiments in conjunction with the accompanying drawings, in which:

[0021] Figure 1 A schematic flowchart of a seismic inversion method according to an exemplary embodiment of the present invention is shown;

[0022] Figure 2 A flowchart illustrating another seismic inversion method according to an exemplary embodiment of the present invention is shown;

[0023] Figure 3 A schematic diagram of a reflection coefficient according to an exemplary embodiment of the present invention is shown;

[0024] Figure 4 A schematic diagram of another reflection coefficient according to an exemplary embodiment of the present invention is shown;

[0025] Figure 5 An exemplary embodiment of the present invention for C is shown. 11 A schematic diagram of contribution analysis;

[0026] Figure 6 An exemplary embodiment of the present invention for C is shown. 13 A schematic diagram of contribution analysis;

[0027] Figure 7 An exemplary embodiment of the present invention for C is shown. 33 A schematic diagram of contribution analysis;

[0028] Figure 8 An exemplary embodiment of the present invention for C is shown.44 A schematic diagram of contribution analysis;

[0029] Figure 9 An exemplary embodiment of the present invention for C is shown. 55 A schematic diagram of contribution analysis;

[0030] Figure 10 A schematic diagram of a contribution analysis for ρ according to an exemplary embodiment of the present invention is shown;

[0031] Figure 11 A schematic diagram of an inversion result according to an exemplary embodiment of the present invention is shown;

[0032] Figure 12 A schematic diagram of another inversion result according to an exemplary embodiment of the present invention is shown;

[0033] Figure 13 An anisotropic parameter δ is shown according to an exemplary embodiment of the present invention. (v) A schematic diagram;

[0034] Figure 14 Another anisotropy parameter δ according to an exemplary embodiment of the present invention is shown. (v) A schematic diagram;

[0035] Figure 15 An anisotropy parameter ε is shown according to an exemplary embodiment of the present invention. (v) A schematic diagram;

[0036] Figure 16 Another anisotropy parameter ε according to an exemplary embodiment of the present invention is shown. (v) A schematic diagram;

[0037] Figure 17 An anisotropic parameter γ is shown according to an exemplary embodiment of the present invention. (v) A schematic diagram;

[0038] Figure 18 Another anisotropy parameter γ according to an exemplary embodiment of the present invention is shown. (v) A schematic diagram. Detailed Implementation

[0039] Embodiments of the present invention will now be described in more detail with reference to the accompanying drawings. While some embodiments of the invention are shown in the drawings, it should be understood that the invention can be implemented in various forms and should not be construed as limited to the embodiments set forth herein. Rather, these embodiments are provided to provide a more thorough and complete understanding of the invention. It should be understood that the accompanying drawings and embodiments are for illustrative purposes only and are not intended to limit the scope of protection of the invention.

[0040] It should be understood that the various steps described in the method embodiments of the present invention may be performed in different orders and / or in parallel. Furthermore, the method embodiments may include additional steps and / or omit the steps shown. The scope of the present invention is not limited in this respect.

[0041] The term "comprising" and its variations as used herein are open-ended, meaning "including but not limited to". The term "based on" means "at least partially based on". The term "one embodiment" means "at least one embodiment"; the term "another embodiment" means "at least one additional embodiment"; the term "some embodiments" means "at least some embodiments". Definitions of other terms will be given in the following description. It should be noted that the concepts of "first", "second", etc., mentioned in this invention are used only to distinguish different devices, modules, or units, and are not intended to limit the order of functions performed by these devices, modules, or units or their interdependencies.

[0042] It should be noted that the terms "a" and "a plurality of" used in this invention are illustrative rather than restrictive. Those skilled in the art should understand that, unless otherwise expressly indicated in the context, they should be understood as "one or more".

[0043] The names of the messages or information exchanged between the multiple devices in the embodiments of the present invention are for illustrative purposes only and are not intended to limit the scope of these messages or information.

[0044] It should be noted that the execution subject of the seismic inversion method provided in this embodiment of the invention can be one or more electronic devices, and this invention does not limit this; wherein, the electronic device can be a terminal (i.e., a client) or a server. Therefore, when the execution subject includes multiple electronic devices, and among these multiple electronic devices are at least one terminal and at least one server, the seismic inversion method provided in this embodiment of the invention can be jointly executed by the terminal and the server. Accordingly, the terminal mentioned herein may include, but is not limited to: smartphones, laptops, desktop computers, intelligent voice interaction devices, etc. The server mentioned herein can be an independent physical server, a server cluster or distributed system composed of multiple physical servers, or a cloud server providing basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communication, middleware services, domain name services, security services, CDN (Content Delivery Network), and big data and artificial intelligence platforms, etc.

[0045] Based on the above description, this embodiment of the invention proposes a seismic inversion method, which can be executed by the aforementioned electronic device (terminal or server); or, the seismic inversion method can be executed jointly by a terminal and a server. For ease of explanation, the following description will use the execution of the seismic inversion method by an electronic device as an example; such as Figure 1 As shown, the seismic inversion method may include the following steps S101-S105:

[0046] S101, acquire initial inversion data, and acquire first azimuth seismic data and second azimuth seismic data.

[0047] Optionally, the initial inversion data can be initial inversion data for any region; correspondingly, the first azimuth seismic data and the second azimuth seismic data can be first azimuth seismic data and second azimuth seismic data for any region. Optionally, the initial inversion data may include the initial values ​​of each of the multiple reference parameters in each sampling medium layer in at least one sampling medium layer, or it may include the initial inversion parameter data (i.e., initial values) of the parameter to be inverted, etc.; this embodiment of the invention does not limit this. For ease of explanation, the following description will use the example that the initial inversion data may include the initial values ​​of each reference parameter in each sampling medium layer. Optionally, the multiple reference parameters may include, but are not limited to: shear wave velocity (V S ), longitudinal wave velocity (V P ), density (ρ) and anisotropy parameter (δ) (v) ε (v) γ (v) )wait.

[0048] Optionally, the initial inversion data, the first azimuth seismic data, and the second azimuth seismic data can be data from the well logging model (i.e., data generated during model testing), or data from actual seismic data, etc.; this embodiment of the invention does not limit this. Optionally, the methods for obtaining the initial inversion data, the first azimuth seismic data, and the second azimuth seismic data can include, but are not limited to, the following:

[0049] The first acquisition method: The electronic device can store the initial inversion data, first azimuth seismic data and second azimuth seismic data of each region in the storage space of multiple regions. In this case, the initial inversion data, first azimuth seismic data and second azimuth seismic data of any region can be obtained from its own storage space.

[0050] The second acquisition method: Electronic equipment can acquire well logging data in any area, and use this well logging data to obtain initial inversion data, first azimuth seismic data and second azimuth seismic data.

[0051] The third acquisition method: The electronic device can display a data input interface. In this case, the user can perform a data input operation on the data input interface to input the initial inversion data, the first azimuth angle seismic data, and the second azimuth angle seismic data. In this case, the electronic device can respond to the data input operation and acquire the initial inversion data, the first azimuth angle seismic data, and the second azimuth angle seismic data, etc.

[0052] S102, construct at least one first azimuth angle forward operator and at least one second azimuth angle forward operator.

[0053] It should be noted that, based on the stress-strain relationship (i.e., the stress-strain relationship in rock), this invention establishes precise elastic impedance tensors for P-waves (longitudinal waves), SV waves (waves where particle vibration occurs in a plane perpendicular to the wave propagation plane), and SH waves (waves where particle vibration occurs in a plane parallel to the wave propagation plane) in HTI media. The elastic impedance tensor, defined by the stress and velocity fields, can completely describe the kinematic characteristics of the rock and accurately represent the reflection and transmission coefficients at arbitrary incident angles and elastic contrast interfaces. In a coordinate system with the z-axis pointing vertically downwards, the stress-strain relationship can be expressed as shown in Equation 1.1:

[0054] Formula 1.1

[0055] Where σ and u represent stress and displacement (σ xz σ represents the stress in the x and z directions. yz σ represents the stress in the y-direction and the z-direction. zz U represents the stress in the z-direction and the z-direction. x u represents the displacement in the x-direction. y u represents the displacement in the y-direction. z (representing displacement in the z-direction), C represents the stiffness matrix of the medium, which can be defined by the stress-strain relationship. ab The elements in the stiffness matrix (also called stiffness coefficients, where a and b are positive integers) can be represented. A stiffness matrix can include the values ​​of each stiffness coefficient. In this embodiment of the invention, the HTI medium requires five independent elastic parameters to characterize it, and its stiffness coefficients can be expressed as shown in Formula 1.2:

[0056] Equation 1.2

[0057] Among them, V S V represents the transverse wave velocity. P ρ represents the longitudinal wave velocity, δ represents the density, and ρ represents the longitudinal wave velocity. (v) ε (v) γ (v)This represents the anisotropy parameter. Correspondingly, based on the displacement equation for seismic wave propagation, Equation 1.1 can be written in matrix-vector form, as shown in Equation 1.3:

[0058] Equation 1.3

[0059] Where i is the imaginary unit, w is the angular frequency, and τ is the phase function, whose gradient is the normal direction of the wavefront. The velocity field can also be expressed in terms of displacement as shown in Equation 1.4:

[0060] Equation 1.4

[0061] Where v is the propagation velocity of the seismic wave in the medium, and t represents time. Based on this, Equation 1.3 represents the stress field related to displacement, and Equation 1.4 represents the characteristics related to the velocity field. The stress field, velocity field, and elastic impedance tensor are related as shown in Equation 1.5:

[0062] Formula 1.5

[0063] Where E represents the elastic impedance tensor of the medium. Based on this, by combining equations 1.3-1.5, the specific expression for the elastic impedance tensor of the HTI medium can be obtained as shown in equation 1.6:

[0064] Equation 1.6

[0065] In this context, the subscript m (m=P,SV,SH) can represent P-wave, SV-wave, or SH-wave; "+" and "-" indicate downward and upward propagation of the seismic wave, i.e., downward and upward waves; p x and p y Let x and y represent the horizontal slowness along the x and y axes, respectively, and q represent the vertical slowness, as shown in Formula 1.7:

[0066] Equation 1.7

[0067] In summary, Equation 1.6 provides the accurate elastic impedance tensors for P-waves, SV-waves, and SH-waves in the HTI medium, which can be expressed using stiffness coefficients and slowness. Based on this, the elastic impedance tensor, defined by the stress and velocity fields, can completely describe the kinematic characteristics of the rock and accurately represent the reflection and transmission coefficients at arbitrary incident angles and elastic contrast interfaces. Furthermore, based on the accurate elastic impedance tensor of the HTI medium, this embodiment derives the accurate and weakly contrast approximate expressions for the reflection coefficient characterized by the stiffness matrix, constructs the propagation matrix for seismic wave propagation in the HTI medium, simplifies the accurate reflection coefficient equation using the relationship between the up and down propagation matrices, and derives an approximate expression for the HTI medium reflection coefficient under the weak contrast assumption, improving the accuracy of the reflection coefficient at large angles. Correspondingly, the polarization vector of the HTI medium can be expressed as shown in Equation 1.8:

[0068] Formula 1.8

[0069] Where U is the polarization vector, G P and G SV The coefficients of the polarization vector can be specifically expressed as shown in Formula 1.9:

[0070] Equation 1.9

[0071] Wherein, the subscript P(SV) indicates that the subscript is P (representing P-wave) or SV (representing SV-wave). Correspondingly, the propagation matrix is ​​an important function describing the propagation characteristics of seismic waves in the medium, and can be expressed by the elastic impedance tensor and polarization vector as shown in Equation 1.10:

[0072] Formula 1.10

[0073] Where K is the HTI medium propagation matrix. In this embodiment of the invention, the relationship between the uplink and downlink elastic impedance tensors and polarization vectors of the HTI medium can be considered as shown in Equation 1.11:

[0074] Equation 1.11

[0075] Where matrix F is a diagonal matrix. Substituting Equation 1.11 into Equation 1.10, the relationship between the uplink and downlink propagation matrices is shown in Equation 1.12:

[0076] Equation 1.12

[0077] Based on this, the propagation matrix and polarization vector can be used to characterize the accurate reflection coefficient equation of the HTI medium, as shown in Equation 1.13:

[0078] Equation 1.13

[0079] In this context, superscripts (1) and (2) represent the upper and lower media, and R is the reflection coefficient composed of qP wave (i.e., P wave), qSV wave (i.e., SV wave), and qSH wave (i.e., SH wave), which can be specifically expressed as shown in Formula 1.14:

[0080] Equation 1.14

[0081] Among them, AB is used to indicate AB waves (A and B are P, SV, or SH respectively), that is, the incident wave is A wave and the reflected wave is B wave. For example, PP is used to indicate PP waves (that is, the incident wave is P wave and the reflected wave is also P wave), PSV is used to indicate PSV waves (that is, the incident wave is P wave and the reflected wave is SV wave), and so on; R AB Let be the reflection coefficient of the AB wave. Based on this, substituting Equation 1.12 into Equation 1.13, the exact reflection coefficient equation can be simplified to Equation 1.15:

[0082] Equation 1.15

[0083] The matrices S and Q can be represented by the propagation matrix and polarization vector as shown in Equation 1.16:

[0084] Equation 1.16

[0085] Here, the superscript n (n=1,2) represents the upper or lower medium. It is assumed that the elastic and anisotropic parameters of the upper and lower media vary relatively little, i.e., Q=(Q (1) +Q (2) ) / 2, S=(S (1) +S (2) ) / 2, ∆Q=Q (2) -Q (1) ∆S=S (2) -S (1) If this holds true, then Equation 1.15 can be further simplified to Equation 1.17:

[0086] Equation 1.17

[0087] Based on this, according to the elastic dynamics equations of seismic wave propagation, the matrices S and Q have a relationship Q T S=S T Q = Λ, where Λ is the scale factor. Therefore, the reflection coefficient of PP waves in HTI medium can be approximately expressed as shown in Equation 1.18:

[0088] Equation 1.18

[0089] Among them, Q P S P Λ PAs shown in Formula 1.19:

[0090] Equation 1.19

[0091] Based on this, substituting Formula 1.19 into Formula 1.18, we can obtain the display form of the PP wave reflection coefficient of the HTI medium as shown in Formula 1.20:

[0092] Equation 1.20

[0093] Among them, v 11 v 13 v 33 v 44 and v 55 As shown in Formula 1.21:

[0094] Equation 1.21

[0095] The first coefficient term Λ1 and the second coefficient term Λ2 can be represented as shown in Formula 1.22:

[0096] Equation 1.22

[0097] Among them, the above Γ 11 =C 11 / C 33 , Γ 13 =C 13 / C 33 , Γ 44 =C 44 / C 33， Γ 55 =C 55 / C 33 , Γ ρ =C ρ / C 33 It should be understood that Equations 1.20 and 1.21 are weakly contrasted approximations of the reflection coefficient of PP waves in HTI media. When the anisotropy parameter is 0, the derived reflection coefficient can be equivalent to the Aki approximation equation in isotropic media.

[0098] Furthermore, this embodiment of the invention establishes forward modeling models for the crack plane (i.e., the first azimuth angle 90°) and the symmetry axis plane (i.e., the second azimuth angle 0°), constructs objective functionals (i.e., the first objective inversion function and the second objective inversion function) that consider the mapping relationship between the stiffness matrix and the anisotropic parameters; and analyzes the variation law of the contribution of each parameter to the reflection coefficient with azimuth, determines that the crack plane and the symmetry axis plane are the dominant azimuths for seismic inversion, and introduces the mapping relationship between the stiffness matrix and the anisotropic parameters, overcoming the problem of ill-posedness in calculating anisotropic parameters using the stiffness matrix.

[0099] It should be noted that, according to the derived equations for the reflection coefficient of PP waves in HTI media (Equations 1.20 and 1.21), when the azimuth (i.e., the azimuth angle) is 90°, the slowness p x If the coefficient is 0, then the coefficient term is as shown in formula 1.23:

[0100] Equation 1.23

[0101] Based on this, parameter C 11 C 13 and C 55 The reflection coefficient does not contribute within the crack plane, meaning that only C needs to be considered within the crack plane. 33 C 44 By using the three parameters ρ, the number of parameters to be inverted is reduced from 6 to 3. Meanwhile, C... 11 C 13 and C 55 The contribution to the reflection coefficient is greatest in the plane of symmetry. Therefore, embodiments of the present invention can use a step-by-step inversion method to predict the stiffness matrix and density: the first step is to predict C in the crack plane. 33 C 44 And ρ; the second step is to predict C in the plane of symmetry. 11 C 13 and C 55 .

[0102] Furthermore, the perturbation terms of each parameter in the reflection coefficient can be written in the form of logarithmic differences, as shown in Equation 1.24:

[0103] Equation 1.24

[0104] Accordingly, based on the reflection coefficient equation, the forward model in the crack plane is established as shown in Equation 1.25:

[0105] Equation 1.25

[0106] Among them, superscript 1 = 90° is the first azimuth angle. This represents seismic data with an azimuth of 90° (also known as first azimuth seismic data), consisting of h incident angles θ1, θ2…θ under the first azimuth. h The Tao set is composed of; X1=[X 33 ,X 44 ,X ρ ] T It can be the vector form of the first parameter to be inverted, determined by parameter C. 33 C 44The logarithm of ρ is used; optionally, the number of first forward operators in at least one first forward operator can be one; optionally, the first forward operator can be expressed as shown in Formula 1.26:

[0107] Equation 1.26

[0108] Wherein, W is the seismic wavelet kernel matrix (also referred to here as the first seismic wavelet kernel matrix), which can be obtained by time shifting the seismic wavelet, such as a wavelet extracted from seismic data at the first azimuth angle; L can be a difference matrix, which can be set according to experience or actual needs, and this embodiment of the invention does not limit this; all parameters before WL can be referred to as the first azimuth angle forward modeling operator perturbation data, that is, the first azimuth angle forward modeling operator perturbation data can include each first perturbation parameter v. 33 v 44 and v ρ The disturbance values ​​at various incident angles and the first azimuth angle; optionally, the disturbance value of a disturbance parameter at an incident angle and an azimuth angle may be affected only by the corresponding incident angle and azimuth angle, that is, other parameters can be canceled out by the above formula to determine the disturbance value of the disturbance parameter, or the value of a disturbance parameter at an incident angle and an azimuth angle may be set according to experience or actual needs.

[0109] Based on this, when constructing at least one first azimuth forward modeling operator, the electronic device can determine the first azimuth forward modeling operator perturbation data and the first seismic wavelet kernel matrix; and construct at least one first azimuth forward modeling operator based on the first azimuth forward modeling operator perturbation data and the first seismic wavelet kernel matrix. Optionally, when determining the first azimuth forward modeling operator perturbation data, h incident angles under the first azimuth can be determined, and the first azimuth forward modeling operator perturbation data can be determined based on the first azimuth and the h incident angles under the first azimuth; the first azimuth forward modeling operator perturbation data includes each first parameter (i.e., X) in the first inversion parameters. 33 X 44 and X ρ The corresponding perturbation parameters (i.e., each of the first perturbation parameters, such as X) 33 X 44 and X ρ respectively with v 33 v 44 and v ρ Correspondingly, the disturbance values ​​at the first azimuth angle and at h incident angles at the first azimuth angle, where h is a positive integer. Therefore, the electronic device can use Formula 1.26 to construct at least one forward operator for the first azimuth angle.

[0110] Similarly, the forward model in the plane of symmetry is shown in Equation 1.27:

[0111] Equation 1.27

[0112] Among them, superscript 0 = 0° is the second azimuth angle. Seismic data with an azimuth of 0° (also known as second azimuth seismic data); X2=[X 11 ,X 13 ,X 55 ] T It can be the vector form of the second parameter to be inverted, determined by parameter C. 11 C 13 and C 55 The logarithmic composition; optionally, the number of second forward operators in at least one second forward operator can be 2, that is, including second forward operators. and Optionally, the second forward operator can be as shown in Equation 1.28:

[0113] Equation 1.28

[0114] Here, W can also be called the second seismic wavelet kernel matrix, which can be a wavelet extracted from the second azimuth seismic data, or a wavelet extracted from the first azimuth seismic data, etc.; optionally, the second seismic wavelet kernel matrix and the first seismic wavelet kernel matrix can be the same or different, and this embodiment of the invention does not limit this; optionally, the second azimuth forward modeling perturbation data can include various second perturbation parameters v. 11 v 13 and v 55 The disturbance values ​​at each incident angle and the second azimuth angle; optionally, the h incident angles at the first azimuth angle and the h incident angles at the second azimuth angle can be the same or different, and the embodiments of the present invention do not limit this; for ease of explanation, the following description will take the example that the h incident angles at the first azimuth angle and the h incident angles at the second azimuth angle are the same.

[0115] Based on this, when constructing at least one second azimuth forward modeling operator, the electronic device can determine the second azimuth forward modeling operator perturbation data and the second seismic wavelet kernel matrix; and construct at least one second azimuth forward modeling operator based on the second azimuth forward modeling operator perturbation data and the second seismic wavelet kernel matrix. Optionally, when determining the second azimuth forward modeling operator perturbation data, h incident angles under the second azimuth are determined, and the second azimuth forward modeling operator perturbation data is determined based on the second azimuth and the h incident angles under the second azimuth; the second azimuth forward modeling operator perturbation data includes the perturbation values ​​of the perturbation parameters corresponding to each first parameter under the second azimuth and the h incident angles under the second azimuth, and each second parameter (i.e., X) in the second parameters to be inverted.11 X 13 and X 55 The corresponding perturbation parameters (i.e., each of the second perturbation parameters, such as X) 11 X 13 and X 55 respectively with v 11 v 13 and v 55 Correspondingly, the disturbance values ​​at the second azimuth angle and the h incident angles at the second azimuth angle. Therefore, the electronic device can use Formula 1.28 to construct at least one second azimuth angle forward operator.

[0116] S103, construct a first objective inversion function based on at least one first azimuth forward modeling operator and first azimuth seismic data; and construct a second objective inversion function based on at least one second azimuth forward modeling operator and second azimuth seismic data.

[0117] S104. Based on the first objective inversion function, the second objective inversion function, and the initial inversion data, determine the objective inversion result of the parameters to be inverted.

[0118] Optionally, when there are multiple sampling medium layers in at least one sampling medium layer, the target inversion result of the parameter to be inverted (i.e., X) may include the target inversion results of each parameter in the parameter to be inverted (i.e., each first parameter in the first parameter to be inverted and each second parameter in the second parameter to be inverted) under each sampling medium layer.

[0119] S105, Calculate the target inversion result of anisotropy parameters based on the target inversion result of the parameters to be inverted.

[0120] Optionally, the target inversion results of the anisotropic parameters may include the target inversion results of the anisotropic parameters under each sampling medium layer.

[0121] This invention, after acquiring initial inversion data and first and second azimuth seismic data, constructs at least one first azimuth forward modeling operator and at least one second azimuth forward modeling operator. Then, based on at least one first azimuth forward modeling operator and the first azimuth seismic data, a first target inversion function is constructed; and based on at least one second azimuth forward modeling operator and the second azimuth seismic data, a second target inversion function is constructed. Further, based on the first and second target inversion functions and the initial inversion data, the target inversion results for the parameters to be inverted are determined; and based on the target inversion results for the parameters to be inverted, the target inversion results for anisotropy parameters are calculated. Therefore, this invention can achieve seismic inversion through the first and second target inversion functions to predict the anisotropy parameters of HTI media, thereby obtaining highly accurate target inversion results for anisotropy parameters and effectively improving the accuracy of seismic inversion.

[0122] Based on the above description, this embodiment of the invention also proposes a more specific seismic inversion method. Accordingly, this seismic inversion method can be executed by the aforementioned electronic device (terminal or server); or, the seismic inversion method can be executed jointly by the terminal and the server. For ease of explanation, the following description will use the execution of this seismic inversion method by an electronic device as an example; please refer to [link to relevant documentation]. Figure 2 The seismic inversion method may include the following steps S201-S207:

[0123] S201, acquire initial inversion data, and acquire first azimuth seismic data and second azimuth seismic data.

[0124] S202, construct at least one first azimuth angle forward operator and at least one second azimuth angle forward operator.

[0125] S203, construct a first objective inversion function based on at least one first azimuth forward modeling operator and first azimuth seismic data; and construct a second objective inversion function based on at least one second azimuth forward modeling operator and second azimuth seismic data.

[0126] In this embodiment of the invention, Equations 1.25 and 1.27 show the forward model in the crack plane and symmetry axis plane of the HTI medium. This forward model can predict the stiffness matrix using a step-by-step inversion method. However, calculating anisotropy using the stiffness matrix is ​​ill-posed; that is, small perturbations in the stiffness matrix can severely affect the prediction of anisotropy parameters. Therefore, this invention introduces a mapping relationship between the stiffness matrix and anisotropy parameters in the inversion, which can be expressed as a function as shown in Equation 2.1:

[0127] Equation 2.1

[0128] Here, the function g can represent the mapping relationship between the stiffness matrix and the anisotropic parameters (i.e., the anisotropic parameters of the HTI medium). Optionally, the parameters to be inverted may include a first parameter to be inverted and a second parameter to be inverted; based on this, when constructing the first target inversion function based on at least one first azimuth forward modeling operator and first azimuth seismic data, the electronic device can construct the first target inversion function based on at least one first azimuth forward modeling operator, first azimuth seismic data, and the second parameter to be inverted. Optionally, the electronic device can use Formula 2.2 to construct the first target inversion function:

[0129] Equation 2.2

[0130] Where k represents the iteration number (i.e., in the k-th iteration, the inversion parameter data of the second parameter to be inverted under the k-th iteration is used as the value of the second parameter to be inverted in the first objective inversion function, so as to predict the inversion parameter data of the first parameter to be inverted under the (k+1)-th iteration), λ1 is the first regularization parameter, which is used to control the smoothness of the stiffness matrix and density inversion results; optionally, the first regularization parameter and parameter α can be set according to experience or actual needs, and the embodiments of the present invention do not limit this; optionally, the first objective inversion function can represent the stiffness matrix prediction in the crack plane.

[0131] Accordingly, when constructing the second objective inversion function based on at least one second azimuth forward modeling operator and second azimuth seismic data, the electronic device can construct the second objective inversion function based on at least one second azimuth forward modeling operator, the second azimuth seismic data, and the first parameter to be inverted. Optionally, the electronic device can use Formula 2.3 to construct the second objective inversion function:

[0132] Equation 2.3

[0133] Wherein, λ2 is the second regularization parameter; optionally, the second regularization parameter can be set according to experience or actual needs, and this embodiment of the invention does not limit this; optionally, the second objective inversion function can represent the prediction of the stiffness matrix in the axial plane. It can be seen that this embodiment of the invention can use the inversion parameter data of the first parameter to be inverted in the (k+1)th iteration as the value of the first parameter to be inverted in the second objective inversion function, to predict the inversion parameter data of the second parameter to be inverted in the (k+1)th iteration. Based on this, the inversion parameter data of the first parameter to be inverted in the (k+1)th iteration can be inverted using Formula 2.2, then the inversion parameter data of the first parameter to be inverted in the (k+1)th iteration can be substituted into Formula 2.3 to invert the inversion parameter data of the second parameter to be inverted in the (k+1)th iteration, and then the inversion parameter data of the second parameter to be inverted in the (k+1)th iteration can be substituted into Formula 2.2 to invert the inversion parameter data of the first parameter to be inverted in the (k+2)th iteration, so as to repeatedly alternately invert the first parameter to be inverted and the second parameter to be inverted until convergence.

[0134] S204. Based on the first objective inversion function, the second objective inversion function, and the initial inversion data, determine the objective inversion result of the parameters to be inverted.

[0135] In this embodiment of the invention, a nonlinear inversion optimization algorithm for HTI media based on azimuth amplitude differences is developed. Due to the introduction of a mapping function between the stiffness matrix and anisotropic parameters, the target functional exhibits nonlinear characteristics. In this embodiment, Newton's method can be used to iteratively solve the target functional. According to Newton's method, the update rules for Equations 2.2 and 2.3 are as shown in Equation 2.4:

[0136] Equation 2.4

[0137] The first step length μ1 and the second step length μ2 can be set empirically or according to actual needs; this embodiment of the invention does not limit this. That is, the values ​​of the first step length and the second step length in the k-th iteration can be set empirically or according to actual needs. Optionally, the values ​​of the first step length and the second step length in different iterations can be the same or different; this embodiment of the invention does not limit this. and This can represent first-order and second-order partial derivative operators. Optionally, the first-order and second-order partial derivatives of the first objective inversion function (i.e., the partial derivatives of the first objective inversion function with respect to the first inversion parameter, also known as the partial derivatives of the first objective inversion function with respect to the first inversion parameter) can be expressed as shown in Equation 2.5:

[0138] Equation 2.5

[0139] Where J1 is the Jacobian matrix of function g with respect to vector X1, and H1 is the Hessian matrix of function g with respect to vector X1. Correspondingly, the first and second partial derivatives of the second objective inversion function (i.e., the partial derivatives of the second objective inversion function with respect to the second inversion parameter, also known as the partial derivatives of the second objective inversion function with respect to the second inversion parameter) can be expressed as shown in Equation 2.6:

[0140] Equation 2.6

[0141] Where J2 is the Jacobian matrix of function g with respect to vector X2, and H2 is the Hessian matrix of function g with respect to vector X2. Optionally, J1, J2, H1, and H2 can be as shown in Equation 2.7:

[0142] Equation 2.7

[0143] Based on this, the Jacobian matrix and Hessian matrix of the function can be obtained according to Formula 2.7, and then the first and second derivatives of the objective functional (i.e., the first objective inversion function and the second objective inversion function) can be calculated. The first and second parameters to be inverted are iteratively inverted based on the update rule.

[0144] In this embodiment of the invention, the electronic device can determine the inversion parameter data of the first parameter to be inverted in the 0th iteration (i.e., the initial inversion parameter data of the first parameter to be inverted) and the inversion parameter data of the second parameter to be inverted in the 0th iteration (i.e., the initial inversion parameter data of the second parameter to be inverted) based on the initial inversion data. For example, when the initial inversion data includes the initial values ​​of each reference parameter under each sampling medium layer, an initial stiffness matrix can be determined based on the initial inversion data (e.g., determined using Formula 1.2, which may include the initial stiffness values ​​of each stiffness coefficient under each sampling medium layer). Based on the initial stiffness matrix, the initial inversion parameter data of the first parameter to be inverted and the initial inversion parameter data of the second parameter to be inverted can be determined (e.g., by performing a logarithmic operation on the stiffness values ​​of the stiffness coefficients). A stiffness matrix may include the stiffness data of each stiffness coefficient among multiple stiffness coefficients, and the stiffness data of a stiffness coefficient may include the stiffness value of the corresponding stiffness coefficient under each sampling medium layer. Based on this, an inversion parameter data may include the inversion parameter values ​​of the first or second parameter to be inverted under each sampling medium layer in one iteration.

[0145] Accordingly, the electronic device can determine the inversion parameter data of the first inversion parameter and the second inversion parameter in the kth iteration based on the first target inversion function, the second target inversion function, and the initial inversion data. In other words, the initial inversion parameter data of the first inversion parameter and the second inversion parameter can be determined based on the initial inversion data. Therefore, based on the first target inversion function, the second target inversion function, the initial inversion parameter data of the first inversion parameter, and the initial inversion parameter data of the second inversion parameter, the inversion parameter data of the first inversion parameter and the second inversion parameter in the kth iteration can be determined, where k is a positive integer.

[0146] Based on this, the partial derivative of the first target inversion function with respect to the first target inversion parameter (i.e., the partial derivative of the first target inversion function with respect to the first target inversion parameter) can be calculated using the inversion parameter data of the first target inversion parameter in the kth iteration and the inversion parameter data of the second target inversion parameter in the kth iteration. For example, the partial derivative of the first target inversion function with respect to the first target inversion parameter (which may include the first-order partial derivative and the second-order partial derivative) can be calculated using Formula 2.5. Based on the partial derivative of the first target inversion function with respect to the first target inversion parameter, the inversion parameter data of the first target inversion parameter in the (k+1)th iteration can be determined. For example, the first target inversion parameter can be updated using Formula 2.4 to obtain the inversion parameter data of the first target inversion parameter in the (k+1)th iteration.

[0147] Accordingly, based on the inversion parameter data of the second inversion parameter in the kth iteration and the inversion parameter data of the first inversion parameter in the (k+1)th iteration, the partial derivative of the second objective inversion function with respect to the second inversion parameter can be calculated (i.e., the partial derivative of the second objective inversion function with respect to the second inversion parameter). For example, formula 2.6 can be used to calculate the partial derivative of the second objective inversion function with respect to the second inversion parameter (which may include the first-order and second-order partial derivatives). Based on the partial derivative of the second objective inversion function with respect to the second inversion parameter, the inversion parameter data of the second inversion parameter in the (k+1)th iteration can be determined. For example, formula 2.4 can be used to update the second inversion parameter to obtain the inversion parameter data of the second inversion parameter in the (k+1)th iteration.

[0148] Furthermore, if the inversion convergence condition is not met, k is incremented by 1 to make k equal to k+1, and k=k+1 is executed. The process continues to iterate based on the inversion parameter data of the first inversion parameter in the kth iteration and the inversion parameter data of the second inversion parameter in the kth iteration. The partial derivative of the first objective inversion function under the first inversion parameter is calculated, and the first and second inversion parameters are iteratively updated until the inversion convergence condition is met. If the inversion convergence condition is met, the inversion parameter data of the first inversion parameter in the (k+1)th iteration is used as the objective inversion result of the first inversion parameter, and the inversion parameter data of the second inversion parameter in the (k+1)th iteration is used as the objective inversion result of the second inversion parameter, so as to determine the objective inversion result of the inversion parameter.

[0149] Optionally, the electronic device can determine the second parameter update difference based on the inversion parameter data of the second parameter to be inverted in the (k+1)th iteration and the second parameter to be inverted in the kth iteration, and determine whether the inversion convergence condition has been met based on the second parameter update difference; alternatively, the electronic device can use Formula 2.8 to determine whether the inversion convergence condition has been met.

[0150] Equation 2.8

[0151] The update difference threshold η can be set based on experience or actual needs, and this embodiment of the invention does not limit this. Therefore, when the update difference of the second parameter is greater than the update difference threshold, it can be determined that the inversion convergence condition has not been met; when the update difference of the second parameter is less than or equal to the update difference threshold, it can be determined that the inversion convergence condition has been met.

[0152] Optionally, the electronic device may also determine the first parameter update difference based on the inversion parameter data of the first parameter to be inverted in the (k+1)th iteration and the inversion parameter data of the first parameter to be inverted in the kth iteration, and determine whether the inversion convergence condition is met based on the first parameter update difference, etc.; the embodiments of the present invention do not limit this.

[0153] Optionally, the electronic device may determine that the inversion convergence condition has been met when the number of iterations reaches a preset iteration threshold; and / or, determine that the inversion convergence condition has not been met when the number of iterations has not reached the preset iteration threshold, etc.; the embodiments of the present invention do not limit this. The optional preset iteration threshold may be set according to experience or according to actual needs, and the embodiments of the present invention do not limit this.

[0154] S205. Based on the target inversion results of the parameters to be inverted, determine the target stiffness matrix. The target stiffness matrix includes the stiffness data of each stiffness coefficient among multiple stiffness coefficients.

[0155] In this embodiment of the invention, the electronic device can perform an exponential operation on the target inversion result of the inversion parameter to be inverted (i.e., perform an exponential operation on the target inversion result of the inversion parameter to be inverted through the exp function) to obtain the target stiffness matrix, thereby obtaining the stiffness value of each stiffness coefficient under each sampling medium layer (which is the target stiffness value here).

[0156] S206, based on the mapping relationship between each stiffness coefficient and multiple reference parameters, determines the calculation method for each reference parameter among the multiple reference parameters, including anisotropic parameters.

[0157] The function g mentioned above can be used to represent the mapping relationship between each stiffness coefficient and multiple reference parameters. In other words, the terms on the right side of Equation 2.1 can be used to calculate V in sequence. P V S , ρ, δ (v) ε (v) γ (v) This represents the mapping relationship between each stiffness coefficient and multiple reference parameters, thus allowing the calculation method for each reference parameter among the multiple reference parameters to be determined. Based on this, electronic equipment can use Formula 2.1 to determine the calculation method for each reference parameter among the multiple reference parameters.

[0158] S207. Based on the calculation methods corresponding to each reference parameter and the target stiffness matrix, determine the target inversion results of the anisotropic parameters.

[0159] In this embodiment of the invention, the electronic device can calculate the target inversion results of each reference parameter according to the calculation method corresponding to each reference parameter and the target stiffness matrix, thereby obtaining the target inversion results of the anisotropic parameters. Based on this, this embodiment of the invention can obtain the target inversion results of the anisotropic parameters under each sampling medium layer.

[0160] Optionally, the electronic device can also determine the PP wave reflection coefficient term; and based on the PP wave reflection coefficient term and the target stiffness matrix, determine the target PP wave reflection coefficient, etc. Based on this, when the number of sampling medium layers in at least one sampling medium layer is multiple, for any two sampling medium layers (i.e., any two HTI medium layers, one upper medium and one lower medium), Formula 1.21 can be used to determine the corresponding PP wave reflection coefficient term for any two sampling medium layers based on the target inversion results of each reference parameter under any two sampling medium layers. For example, based on the target inversion results of each reference parameter under any two sampling medium layers, the inversion mean of each reference parameter can be determined, and Formula 1.21 can be used to determine the corresponding PP wave reflection coefficient term for any two sampling medium layers based on the inversion mean of each reference parameter, etc. Optionally, the electronic device can use Formula 1.20 to determine the target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix, where the difference between any stiffness coefficient (i.e., ...) The value can be the difference between the target stiffness value of any stiffness coefficient of the lower medium and the target stiffness value of any stiffness coefficient of the upper medium in any two sampled medium layers. The value of any stiffness coefficient substituted into formula 1.20 can be the average value between the target stiffness values ​​of any stiffness coefficient of the lower medium and the target stiffness values ​​of any stiffness coefficient of the upper medium, and so on.

[0161] In this embodiment of the invention, to further verify the feasibility of the anisotropic parameter prediction method proposed in this embodiment, on the one hand, this embodiment tests the accuracy of the proposed HTI medium reflection coefficient approximation on a two-layer medium. Both the upper and lower layers are HTI media, and the differences in parameters are small. By comparing with the accurate value of the reflection coefficient, the effectiveness and accuracy of the proposed reflection coefficient approximation equation are tested, and the contribution of each parameter to the reflection coefficient with azimuth is analyzed to select the advantageous azimuth for stiffness matrix and anisotropic parameter inversion. On the other hand, the feasibility of the proposed inversion method (i.e., seismic inversion method) is tested on models and actual seismic data. In the model test, noise is added to the synthesized seismic record to test the stability and noise resistance of the inversion method. In the actual seismic data test, by comparing with conventional anisotropic seismic inversion methods, the reliability and accuracy of the proposed inversion in anisotropic prediction of fractured reservoirs are studied.

[0162] For example, when testing the accuracy of the approximate equation for the reflection coefficient on a two-layer medium, assume that the stiffness matrix and density of the upper medium are C0 and C1, respectively. 11 =12.19 GPa (gigapascals), C 13 =4.88Gpa, C 33 =13.54Gpa, C 44 =4.97Gpa, C 55=3.98 GPa, ρ=2.64 GPa, and the anisotropy parameter is δ (v) =-0.05、ε (v) =-0.05、γ (v) =-0.1; the stiffness matrix and density of the lower medium are C 11 =15.19 GPa, C 13 =7.50Gpa, C 33 =16.88Gpa, C 44 =6.08Gpa, C 55 =4.25 GPa, ρ=2.70 GPa, anisotropy parameter is δ (v) =-0.1、ε (v) =-0.1、γ (v) =-0.15. Based on this, Figure 3 It shows different directions (i.e., azimuth angles). The reflection coefficient (with different values) varies with the incident angle (i.e., θ), where the black line is the exact reflection coefficient, the short dashed line is the classic Rüger approximation (a weak anisotropic approximation), and the double dashed line (i.e., the long dashed line) is the reflection coefficient proposed in the embodiments of this invention. PP Indicates the reflection coefficient of PP waves; Figure 4 The variation of the reflection coefficient with azimuth at different incident angles is shown. It can be seen that the reflection coefficient proposed in the embodiments of the present invention matches the accurate value better, indicating that the HTI medium reflection coefficient proposed in the present invention is effective. That is, the reflection coefficient calculated by the target inversion result obtained by the inversion method proposed in the embodiments of the present invention is effective, i.e., the accuracy is higher.

[0163] For example, the contribution of each parameter to the reflection coefficient was analyzed on this two-layer model. Figures 5-10 The stiffness matrix C is shown. 11 C 13 C 33 C 44 C 55 The contribution of density to the reflection coefficient at different orientations varies with the incident angle. It can be seen that the stiffness coefficient C... 11 C 13 and C 55 Within the crack plane ( =90°) contributes nothing to the reflection coefficient, stiffness coefficient C 11 C 13 and C 55 In the plane of symmetry ( =0°) contributes significantly to the reflection coefficient, etc.

[0164] For example, the feasibility of the proposed anisotropic parameter prediction method was tested on a model in this embodiment of the invention. The model is based on well logging curves obtained from a reservoir with vertical fractures. It is worth noting that the anisotropic parameters are calculated by a rock physics model, and seismic records are synthesized using an accurate reflection coefficient equation. By adding a certain amount of noise, noise-free and signal-to-noise ratio (SNR) seismic data are generated. Figure 11 These are the elastic and anisotropic parameters predicted from noiseless seismic data; where km / s represents kilometers per second, and g / cm³ represents the elastic and anisotropic parameters predicted from noiseless seismic data. 3 The gray line represents grams per cubic centimeter; the gray line represents the initial model (i.e., the initial inversion data), obtained by low-pass filtering of the real model; the black line represents the real model; the dotted dashed line represents the inversion result of the conventional linear method; and the dashed line represents the target inversion result predicted by the method proposed in this embodiment of the invention. Under noise-free conditions, the target inversion result predicted by the method proposed in this embodiment of the invention matches the real model well. Correspondingly, Figure 12 The inversion result with a signal-to-noise ratio of 8 is slightly lower in terms of signal-to-noise ratio and accuracy of the anisotropic parameters compared to the noiseless inversion result. However, the method proposed in this embodiment is still superior to the conventional linear method, proving that the inversion method proposed in this invention is effective.

[0165] For example, the proposed inversion method can be applied to a specific work area in an embodiment of the present invention. Figure 13 and Figure 14 The anisotropy parameters predicted by the inversion method proposed in this embodiment and the conventional linear method are shown in the profile. High-value areas indicate reservoir locations with vertical fractures. The black curves (i.e., the longitudinal black curves) are the corresponding anisotropy parameter curves (indicating that curves on the right (to the right of the dashed line) indicate higher values, and curves on the left indicate lower values; that is, the further to the right the black curve is, the larger the parameter value it indicates). Compared with the conventional inversion method, the method proposed in this embodiment has advantages in resolution and accuracy. Accordingly, Figure 15 and Figure 16 These are the anisotropy parameters predicted by the inversion method and conventional linear method proposed in this embodiment of the invention, and... Figure 17 and Figure 18 These are the anisotropy parameters predicted by the inversion method proposed in this embodiment of the invention and the conventional linear method. A comparison of the prediction results of the conventional method and the method proposed in this embodiment shows that the inversion method proposed in this embodiment effectively improves the stability and accuracy of anisotropy parameter prediction.

[0166] In this embodiment of the invention, after acquiring initial inversion data, and after acquiring first azimuth seismic data and second azimuth seismic data, at least one first azimuth forward modeling operator and at least one second azimuth forward modeling operator are constructed. Then, based on at least one first azimuth forward modeling operator and the first azimuth seismic data, a first target inversion function is constructed; and based on at least one second azimuth forward modeling operator and the second azimuth seismic data, a second target inversion function is constructed. Based on this, the target inversion result of the parameters to be inverted is determined based on the first target inversion function, the second target inversion function, and the initial inversion data; and based on the target inversion result of the parameters to be inverted, a target stiffness matrix is ​​determined, the target stiffness matrix including the stiffness data of each stiffness coefficient among multiple stiffness coefficients. Further, based on the mapping relationship between each stiffness coefficient and multiple reference parameters, the calculation method corresponding to each of the multiple reference parameters, including anisotropic parameters, is determined; and according to the calculation method corresponding to each reference parameter and the target stiffness matrix, the target inversion result of the anisotropic parameters is determined. It is evident that the embodiments of the present invention can effectively improve the accuracy of the target inversion results of anisotropy parameters, thereby improving the reliability and accuracy of inversion in anisotropy prediction of fractured reservoirs.

[0167] Based on the description of the relevant embodiments of the above-mentioned seismic inversion method, this invention also proposes a seismic inversion device, which can be a computer program (including program code) running in an electronic device; the seismic inversion device can include an acquisition unit and a processing unit. The seismic inversion device can execute... Figure 1 or Figure 2 The earthquake inversion method shown, i.e., the earthquake inversion device can operate the above-mentioned units:

[0168] The acquisition unit is used to acquire initial inversion data, as well as first azimuth seismic data and second azimuth seismic data;

[0169] The processing unit is used to construct at least one first azimuth angle forward operator and at least one second azimuth angle forward operator;

[0170] The processing unit is further configured to construct a first objective inversion function based on the at least one first azimuth forward modeling operator and the first azimuth seismic data; and to construct a second objective inversion function based on the at least one second azimuth forward modeling operator and the second azimuth seismic data;

[0171] The processing unit is further configured to determine the target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function, and the initial inversion data;

[0172] The processing unit is further configured to calculate the target inversion result of the anisotropy parameters based on the target inversion result of the parameters to be inverted.

[0173] In one implementation, when constructing at least one first azimuth forward modeling operator, the processing unit may specifically be used to: determine the first azimuth forward modeling operator perturbation data and determine the first seismic wavelet kernel matrix; and construct at least one first azimuth forward modeling operator based on the first azimuth forward modeling operator perturbation data and the first seismic wavelet kernel matrix.

[0174] When constructing at least one second azimuth forward modeling operator, the processing unit can specifically be used to: determine the second azimuth forward modeling operator perturbation data and determine the second seismic wavelet kernel matrix; and construct at least one second azimuth forward modeling operator based on the second azimuth forward modeling operator perturbation data and the second seismic wavelet kernel matrix.

[0175] In another embodiment, when determining the first azimuth forward modeling operator perturbation data, the processing unit may specifically be used to: determine h incident angles under the first azimuth, and determine the first azimuth forward modeling operator perturbation data based on the first azimuth and the h incident angles under the first azimuth; the first azimuth forward modeling operator perturbation data includes the perturbation values ​​of the perturbation parameters corresponding to each first parameter in the first parameters to be inverted in the first parameters under the first azimuth and the h incident angles under the first azimuth, where h is a positive integer;

[0176] When determining the perturbation data of the second azimuth forward modeling operator, the processing unit can specifically be used to: determine h incident angles under the second azimuth, and determine the perturbation data of the second azimuth forward modeling operator based on the second azimuth and the h incident angles under the second azimuth; the perturbation data of the second azimuth forward modeling operator includes the perturbation values ​​of the perturbation parameters corresponding to each of the first parameters under the second azimuth and the h incident angles under the second azimuth, and the perturbation values ​​of the perturbation parameters corresponding to each of the second parameters in the second parameters to be inverted under the second azimuth and the h incident angles under the second azimuth.

[0177] In another embodiment, the parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted; when the processing unit constructs the first target inversion function based on the at least one first azimuth forward modeling operator and the first azimuth seismic data, it can be specifically used to: construct the first target inversion function based on the at least one first azimuth forward modeling operator, the first azimuth seismic data and the second parameter to be inverted;

[0178] When the processing unit constructs the second target inversion function based on the at least one second azimuth forward modeling operator and the second azimuth seismic data, it can specifically be used to: construct the second target inversion function based on the at least one second azimuth forward modeling operator, the second azimuth seismic data, and the first parameter to be inverted.

[0179] In another embodiment, the parameters to be inverted include a first parameter to be inverted and a second parameter to be inverted. When the processing unit determines the target inversion result of the parameters to be inverted based on the first target inversion function, the second target inversion function, and the initial inversion data, it can specifically be used to: determine the inversion parameter data of the first parameter to be inverted in the k-th iteration and the inversion parameter data of the second parameter to be inverted in the k-th iteration based on the first target inversion function, the second target inversion function, and the initial inversion data; calculate the partial derivative of the first target inversion function with respect to the first parameter to be inverted, where k is a positive integer, based on the partial derivative of the first target inversion function with respect to the first parameter to be inverted; and determine the inversion parameter data of the first parameter to be inverted in the (k+1)-th iteration based on the partial derivative of the first target inversion function with respect to the first parameter to be inverted; and determine the inversion parameter data of the first parameter to be inverted in the (k+1)-th iteration based on the inversion parameter data of the second parameter to be inverted in the k-th iteration and the first parameter to be inverted. In the (k+1)th iteration, the partial derivative of the second target inversion function with respect to the second inversion parameter is calculated based on the inversion parameter data. Based on this partial derivative, the inversion parameter data of the second inversion parameter in the (k+1)th iteration is determined. If the inversion convergence condition is not met, k is incremented by 1 to make k equal to k+1. The process of calculating the partial derivative of the first target inversion function with respect to the first inversion parameter based on the inversion parameter data of the first inversion parameter in the kth iteration and the inversion parameter data of the second inversion parameter in the kth iteration continues iteratively until the inversion convergence condition is met. If the inversion convergence condition is met, the inversion parameter data of the first inversion parameter in the (k+1)th iteration is used as the target inversion result of the first inversion parameter, and the inversion parameter data of the second inversion parameter in the (k+1)th iteration is used as the target inversion result of the second inversion parameter, thereby determining the target inversion result of the parameter to be inverted.

[0180] In another embodiment, when the processing unit calculates the target inversion result of the anisotropic parameter based on the target inversion result of the parameter to be inverted, it may specifically be used to: determine the target stiffness matrix based on the target inversion result of the parameter to be inverted, the target stiffness matrix including the stiffness data of each stiffness coefficient among a plurality of stiffness coefficients; determine the calculation method corresponding to each of the plurality of reference parameters based on the mapping relationship between each stiffness coefficient and a plurality of reference parameters, the plurality of reference parameters including the anisotropic parameter; and determine the target inversion result of the anisotropic parameter according to the calculation method corresponding to each reference parameter and the target stiffness matrix.

[0181] In another embodiment, the processing unit can also be used to: determine the PP wave reflection coefficient term; and determine the target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix.

[0182] According to one embodiment of the present invention, each unit in the seismic inversion device can be individually or entirely merged into one or more other units, or some of the units can be further divided into multiple functionally smaller units. This achieves the same operation without affecting the technical effect of the embodiments of the present invention. The above-mentioned units are based on logical function division. In practical applications, the function of one unit can also be implemented by multiple units, or the function of multiple units can be implemented by one unit. In other embodiments of the present invention, any seismic inversion device may also include other units. In practical applications, these functions can also be implemented with the assistance of other units, and can be implemented collaboratively by multiple units.

[0183] According to another embodiment of the present invention, it is possible to perform operations such as those described above by running on a general-purpose electronic device, such as a computer, which includes processing elements and storage elements such as a central processing unit (CPU), random access memory (RAM), and read-only memory (ROM). Figure 1 or Figure 2 The computer program (including program code) for each step involved in the corresponding method shown herein is used to construct a seismic inversion device and to implement the seismic inversion method of the embodiments of the present invention. The computer program may be stored on, for example, a computer storage medium, loaded onto the aforementioned electronic device via the computer storage medium, and run therein.

[0184] Based on the description of the method and apparatus embodiments above, an exemplary embodiment of the present invention also provides an electronic device, including: at least one processor; and a memory communicatively connected to the at least one processor. The memory stores a computer program executable by the at least one processor, which, when executed by the at least one processor, causes the electronic device to perform the method according to an embodiment of the present invention.

[0185] An exemplary embodiment of the present invention also provides a non-transitory computer-readable storage medium storing a computer program, wherein the computer program, when executed by a computer's processor, is used to cause the computer to perform a method according to an embodiment of the present invention.

[0186] An exemplary embodiment of the present invention also provides a computer program product, including a computer program, wherein, when executed by a computer's processor, the computer program is used to cause the computer to perform a method according to an embodiment of the present invention.

[0187] Furthermore, it should be understood that the above-disclosed embodiments are merely preferred embodiments of the present invention and should not be construed as limiting the scope of the present invention. Therefore, any equivalent variations made in accordance with the claims of the present invention are still within the scope of the present invention.

Claims

1. A seismic inversion method characterized in that, The method comprises: acquiring initial inversion data, and acquiring first azimuthal seismic data and second azimuthal seismic data; the initial inversion data comprises initial values of each of a plurality of reference parameters under each of a plurality of sampling medium layers, the plurality of reference parameters comprising shear wave velocity, longitudinal wave velocity, density, and anisotropy parameter; constructing at least one first azimuthal forward operator and at least one second azimuthal forward operator; based on the at least one first azimuthal forward operator and the first azimuthal seismic data, constructing a first target inversion function; and based on the at least one second azimuthal forward operator and the second azimuthal seismic data, constructing a second target inversion function; based on the first target inversion function, the second target inversion function, and the initial inversion data, determining a target inversion result of a parameter to be inverted; based on the target inversion result of the parameter to be inverted, calculating a target inversion result of the anisotropy parameter, comprising: based on the target inversion result of the parameter to be inverted, determining a target stiffness matrix, the target stiffness matrix comprising stiffness data of each of a plurality of stiffness coefficients; based on a mapping relationship between each of the plurality of stiffness coefficients and a plurality of reference parameters, including the anisotropy parameter, determining a corresponding calculation mode of each of the plurality of reference parameters; and based on the corresponding calculation mode of each of the plurality of reference parameters and the target stiffness matrix, determining the target inversion result of the anisotropy parameter; wherein the target inversion result of the anisotropy parameter comprises a target inversion result of the anisotropy parameter under each of the plurality of sampling medium layers; The method further comprises: determining a PP wave reflection coefficient term, comprising: for any two of the plurality of sampling medium layers, based on the target inversion result of each of the plurality of reference parameters under the any two of the plurality of sampling medium layers, determining an inversion mean value of each of the plurality of reference parameters, and based on the inversion mean value of each of the plurality of reference parameters, determining a corresponding PP wave reflection coefficient term of the any two of the plurality of sampling medium layers; based on the PP wave reflection coefficient term and the target stiffness matrix, determining a target PP wave reflection coefficient, comprising: based on the corresponding PP wave reflection coefficient term of the any two of the plurality of sampling medium layers and the target stiffness matrix, determining a corresponding target PP wave reflection coefficient of the any two of the plurality of sampling medium layers; wherein the target stiffness matrix comprises target stiffness values of each of the plurality of stiffness coefficients under each of the plurality of sampling medium layers, and the corresponding target PP wave reflection coefficient of the any two of the plurality of sampling medium layers is determined based on a mean value between a target stiffness value of any one of the plurality of stiffness coefficients of a lower medium and a target stiffness value of the any one of the plurality of stiffness coefficients of an upper medium.

2. The method of claim 1, wherein, The method further comprises: determining a first azimuthal forward operator perturbation data and a first seismic wavelet kernel matrix; based on the first azimuthal forward operator perturbation data and the first seismic wavelet kernel matrix, constructing at least one first azimuthal forward operator; The constructing at least one second azimuthal angle forward operator comprises: determining second azimuthal angle forward operator perturbation data, and determining a second seismic wavelet kernel matrix; constructing at least one second azimuthal angle forward operator based on the second azimuthal angle forward operator perturbation data and the second seismic wavelet kernel matrix.

3. The method of claim 2, wherein, The determining first azimuthal angle forward operator perturbation data comprises: determining h incident angles under a first azimuthal angle, and determining first azimuthal angle forward operator perturbation data based on the first azimuthal angle and the h incident angles under the first azimuthal angle; the first azimuthal angle forward operator perturbation data comprises perturbation values of perturbation parameters corresponding to each first parameter in first to-be-inverted parameters under the first azimuthal angle and the h incident angles under the first azimuthal angle, and h is a positive integer; The determining second azimuthal angle forward operator perturbation data comprises: determining h incident angles under a second azimuthal angle, and determining second azimuthal angle forward operator perturbation data based on the second azimuthal angle and the h incident angles under the second azimuthal angle; the second azimuthal angle forward operator perturbation data comprises perturbation values of perturbation parameters corresponding to each first parameter under the second azimuthal angle and the h incident angles under the second azimuthal angle, and perturbation values of perturbation parameters corresponding to each second parameter in second to-be-inverted parameters under the second azimuthal angle and the h incident angles under the second azimuthal angle.

4. The method according to any one of claims 1 to 3, characterized in that, The to-be-inverted parameters comprise first to-be-inverted parameters and second to-be-inverted parameters; the constructing a first target inversion function based on the at least one first azimuthal angle forward operator and the first azimuthal angle seismic data comprises: constructing a first target inversion function based on the at least one first azimuthal angle forward operator, the first azimuthal angle seismic data and the second to-be-inverted parameters; The constructing a second target inversion function based on the at least one second azimuthal angle forward operator and the second azimuthal angle seismic data comprises: constructing a second target inversion function based on the at least one second azimuthal angle forward operator, the second azimuthal angle seismic data and the first to-be-inverted parameters.

5. The method according to any one of claims 1 to 3, characterized in that, The to-be-inverted parameters comprise first to-be-inverted parameters and second to-be-inverted parameters; the determining a target inversion result of to-be-inverted parameters based on the first target inversion function, the second target inversion function and the initial inversion data comprises: determining inversion parameter data of the first to-be-inverted parameters under the kth iteration and inversion parameter data of the second to-be-inverted parameters under the kth iteration based on the first target inversion function, the second target inversion function and the initial inversion data, and calculating a partial derivative of the first target inversion function under the first to-be-inverted parameters based on the inversion parameter data of the first to-be-inverted parameters under the kth iteration and the inversion parameter data of the second to-be-inverted parameters under the kth iteration, and determining inversion parameter data of the first to-be-inverted parameters under the k+1th iteration based on the partial derivative of the first target inversion function under the first to-be-inverted parameters, and k is a positive integer; calculating a partial derivative of the second objective inversion function at the second to-be-inverted parameter based on the inversion parameter data of the second to-be-inverted parameter at the kth iteration and the inversion parameter data of the first to-be-inverted parameter at the k+1th iteration; and determining the inversion parameter data of the second to-be-inverted parameter at the k+1th iteration based on the partial derivative of the second objective inversion function at the second to-be-inverted parameter; if the inversion convergence condition is not reached, increasing k by 1 so that k equals k+1, and continuing to iteratively perform the calculation of the partial derivative of the first objective inversion function at the first to-be-inverted parameter based on the inversion parameter data of the first to-be-inverted parameter at the kth iteration and the inversion parameter data of the second to-be-inverted parameter at the kth iteration until the inversion convergence condition is reached; if the inversion convergence condition is reached, taking the inversion parameter data of the first to-be-inverted parameter at the k+1th iteration as the target inversion result of the first to-be-inverted parameter and taking the inversion parameter data of the second to-be-inverted parameter at the k+1th iteration as the target inversion result of the second to-be-inverted parameter to achieve the determination of the target inversion result of the to-be-inverted parameter.

6. A seismic inversion apparatus characterized by, The apparatus comprises: an acquisition unit configured to acquire initial inversion data, and acquire first azimuthal seismic data and second azimuthal seismic data; the initial inversion data comprises initial values of each reference parameter in each sampling medium layer in a plurality of reference parameters, and the plurality of reference parameters comprises a shear wave velocity, a longitudinal wave velocity, a density, and an anisotropy parameter; a processing unit configured to construct at least one first azimuthal forward operator, and construct at least one second azimuthal forward operator; the processing unit is further configured to construct a first objective inversion function based on the at least one first azimuthal forward operator and the first azimuthal seismic data, and construct a second objective inversion function based on the at least one second azimuthal forward operator and the second azimuthal seismic data; the processing unit is further configured to determine a target inversion result of to-be-inverted parameters based on the first objective inversion function, the second objective inversion function, and the initial inversion data; the processing unit is further configured to calculate a target inversion result of the anisotropy parameter based on the target inversion result of the to-be-inverted parameters, including: determining a target stiffness matrix based on the target inversion result of the to-be-inverted parameters, the target stiffness matrix comprising stiffness data of each stiffness coefficient in a plurality of stiffness coefficients; determining a calculation mode corresponding to each reference parameter in a plurality of reference parameters based on a mapping relationship between the each stiffness coefficient and the plurality of reference parameters, the plurality of reference parameters comprising the anisotropy parameter; and determining the target inversion result of the anisotropy parameter according to the calculation mode corresponding to the each reference parameter and the target stiffness matrix; wherein the target inversion result of the anisotropy parameter comprises a target inversion result of the anisotropy parameter in each sampling medium layer. The processing unit is further configured to determine a PP wave reflection coefficient term, including: determining, for any two sampling medium layers in the plurality of sampling medium layers, an inversion mean value of the respective reference parameter based on a target inversion result of the respective reference parameter under the any two sampling medium layers, and determining a PP wave reflection coefficient term corresponding to the any two sampling medium layers based on the inversion mean value of the respective reference parameter; and determining a target PP wave reflection coefficient based on the PP wave reflection coefficient term and the target stiffness matrix, including: determining a target PP wave reflection coefficient corresponding to the any two sampling medium layers based on the PP wave reflection coefficient term corresponding to the any two sampling medium layers and the target stiffness matrix; wherein the target stiffness matrix comprises target stiffness values of the respective stiffness coefficients under the respective sampling medium layers, and the target PP wave reflection coefficient corresponding to the any two sampling medium layers is determined based on a mean value between a target stiffness value of any stiffness coefficient of a lower medium layer and a target stiffness value of the any stiffness coefficient of an upper medium layer in the any two sampling medium layers.

7. An electronic device, comprising: The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program.

8. A non-transitory computer readable storage medium having stored thereon computer instructions, wherein, The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer program. The computer program product comprises a computer readable medium, and the computer readable medium stores the computer

Images