Biaxial stress induced anisotropic medium azimuth seismic inversion method

Abstract

The invention relates to the field of oil-gas exploration and development, in particular to a biaxial stress induced anisotropic medium azimuth seismic inversion method. The method comprises the following steps: S1, parameterizing the rigidity tensor of the biaxial stress induced anisotropic medium; s2, deducing a biaxial stress induced anisotropic medium PP wave reflection coefficient equation; and S3, biaxial stress induced anisotropic medium azimuth seismic inversion is carried out. According to the method, the model parameters of the anisotropic part are inverted by using the amplitude difference of seismic reflection data in different directions, and then the model parameters of the anisotropic part obtained by inversion are used as input to further invert the model parameters of the isotropic part. According to the method, the morbidity of multi-parameter synchronous inversion is effectively reduced, and the stability and precision of inversion are improved.

Description

Technical Field

[0001] The present invention relates to the field of oil and gas exploration and development, and in particular to a method for azimuthal seismic inversion of biaxial stress-induced anisotropic media. Background Art

[0002] Seismic anisotropy reveals the directional dependence of seismic wave velocity and amplitude, and is of great significance for geophysical exploration and development. The physical mechanisms that cause anisotropy are diverse, including preferentially oriented crystals or minerals, fine-scale stratification or arrangement of clay particles, the accumulation of strata with varying elastic properties, directional fractures, and in-situ stress. Like intrinsic and fracture-induced anisotropy, stress-induced anisotropy is also widespread in the subsurface, yet it is often overlooked in seismic anisotropy studies.

[0003] Currently, the PP wave reflection coefficient equation for stress-induced anisotropic media is typically derived based on acoustoelastic theory. This approach treats the rock as a continuous medium and explores the relationship between the inherent stress state and elastic properties from a macroscopic perspective. For example, the technical solution described in the patent document "A Method and System for Azimuthal Seismic Inversion of Stress-Horizontally Anisotropic Media" (Application Publication No. CN115840250A) is as follows: First, a transversely isotropic medium (HTI) model with a horizontal axis of symmetry under linear slip boundary conditions is established, and a small uniaxial stress parallel to the axis of symmetry is applied to the HTI model. Second, the stiffness tensor induced by fractures and stress in the horizontally uniaxially stressed HTI medium is decoupled. Third, an approximate PP wave equation is established that decouples fractures and stress-induced anisotropy parameters. Fourth, azimuthal seismic inversion is performed using the difference in reflection coefficient amplitudes between different azimuths. This method focuses on HTI media subjected to horizontal stress and utilizes solid acoustoelastic theory to study the effects of stress on rock properties.

[0004] However, the stress dependence of rock elastic properties can be explained as the response of microcracks in the rock to an applied stress field. Specifically, when an isotropic rock is subjected to anisotropic stress, microcracks perpendicular to the stress direction preferentially close, inducing elastic anisotropy in the rock. However, the acoustoelastic theory cannot explain the mechanism of stress-induced elastic anisotropy at the microscopic level. Furthermore, the PP wave reflection coefficient of anisotropic media contains many unknown parameters, and the numerical contributions of different parameters to the reflection coefficient vary significantly.

[0005] Therefore, existing multi-parameter simultaneous inversion of anisotropic media often suffers from problems such as low parameter inversion stability and poor noise resistance. Summary of the Invention

[0006] The present invention takes into account the effect of stress on microcrack closure, and based on this, derives a PP wave reflection coefficient equation for coupling stress and rock microstructure changes. It also proposes an azimuthal seismic step-by-step inversion method based on azimuthal amplitude differences under a Bayesian framework to overcome the above-mentioned defects of the prior art.

[0007] To solve the above technical problems, the present invention adopts a technical solution: a method for azimuthal seismic inversion of biaxial stress-induced anisotropic media, which comprises the following steps:

[0008] Step S1, parameterization of the stiffness tensor of the biaxial stress-induced anisotropic medium;

[0009] Step S2, deriving the PP wave reflection coefficient equation of the biaxial stress-induced anisotropic medium;

[0010] Step S3: azimuthal seismic inversion of biaxial stress-induced anisotropic media.

[0011] Preferably, the step S1 includes:

[0012] Considering that the unstressed rock is elastically isotropic and filled with randomly distributed and oriented microcracks, the elastic stiffness tensor of the unstressed rock is expressed as:

[0013] (1),

[0014] In formula (1), and represent the longitudinal wave modulus and shear modulus of unstressed rock respectively;

[0015] The equivalent compliance tensor caused by the closure of microcracks due to horizontal biaxial stress is expressed as:

[0016] (2),

[0017] In formula (2), represents the equivalent flexibility matrix of rock without microcracks, and They represent the normal flexibility and tangential flexibility of microcracks in unstressed rock, respectively. The normalized stress is expressed as and , The tensor represents the typical closure pressure of the crack, and represents the compliance contribution of the initial microcracks, the tensor and Represents stress The flexibility contribution to crack closure, tensor and Represents stress The tensor expressions of the flexibility contribution that causes crack closure are:

[0018] (3),

[0019] (4),

[0020] (5),

[0021] Rearranging formula (2) yields:

[0022] (6),

[0023] In formula (6), represents the flexibility tensor of the unstressed rock, and Represents horizontal stress 、 Compliance contribution that causes microcrack closure;

[0024] Assuming biaxial stress-induced anisotropic media, the equivalent stiffness tensor is derived as:

[0025] (7),

[0026] Combining formulas (1)-(7), the stiffness tensor of biaxial stress-induced anisotropic media is derived as:

[0027] (8),

[0028] The biaxial stress-induced anisotropic medium is approximately equivalent to the orthotropic medium. The stress-induced anisotropic parameter is defined as:

[0029] (9),

[0030] (10),

[0031] The components of the elastic stiffness tensor in formula (8) are expressed as:

[0032] (11),

[0033] In formula (11), , , .

[0034] Preferably, the step S2 includes:

[0035] The PP wave reflection coefficient equation is expressed as:

[0036] (12),

[0037] In formula (12), represents the density of the background medium, is the angle of incidence, represents the perturbation of the stiffness tensor components, The expression is as follows:

[0038] (13),

[0039] In formula (13), represents the observation azimuth, represents the longitudinal wave velocity of the background medium;

[0040] Consider the interface of two biaxial stress-induced anisotropic media with different properties, and combine formula (11) to derive the perturbation terms of each component of the stiffness tensor, ignoring the , , , , , , , , , , , Proportional terms, we get:

[0041] (14),

[0042] In formula (14), △ represents the difference in the properties of the media above and below the reflection interface;

[0043] Combining formulas (12)-(14), the PP wave reflection coefficient equation of biaxial stress-induced anisotropic media is derived as:

[0044] (15),

[0045] In formula (15),

[0046] (16),

[0047] (17),

[0048] (18),

[0049] In formulas (16)-(18), , , ,

[0050] , ,

[0051] ,

[0052] .

[0053] Preferably, step S3 includes:

[0054] Step S31: According to formula (15), the azimuth difference reflection coefficient is expressed as:

[0055] (19),

[0056] In formula (19), and represent two different observation azimuths, and

[0057] (20),

[0058] Considering the case of M reflection interfaces and L incident angles, based on the convolution model, the azimuth difference of seismic data is expressed in matrix form as follows:

[0059] (twenty one),

[0060] In formula (21),

[0061] ,

[0062] The incident angle is And the observation azimuth is The earthquake trace;

[0063] ,

[0064] ,

[0065] ,

[0066] , , ,

[0067] , ,

[0068] Where W represents the wavelet matrix, represents the Mth reflection interface, and Respectively represent matrix transposition and generate an M×M diagonal matrix;

[0069] In the Bayesian theory framework, the posterior probability distribution of the model parameters is Proportional to the likelihood function of seismic data and the prior probability distribution of the model parameters :

[0070] (twenty two),

[0071] In formula (22), Assuming it obeys the Cauchy distribution, it can be expressed as:

[0072] (twenty three),

[0073] In formula (23), represents the variance of the model parameters;

[0074] Assuming that the likelihood function follows a Gaussian distribution, it can be expressed as:

[0075] (twenty four),

[0076] In formula (24), represents the variance of noise in seismic data;

[0077] Combining formulas (22)-(24), the target functional of step S31 is constructed as:

[0078] (25),

[0079] In formula (25), represents the regularization coefficient of the model parameters, Represents the low-frequency components of the model parameters;

[0080] The iterative weighted least squares algorithm is used to optimize and solve formula (25). After multiple iterations, the unknown model parameters are obtained;

[0081] Step S32: The model parameters of the anisotropic part obtained by inversion in step S31 are , as input, invert the model parameters of the isotropic part;

[0082] The earthquake forward modeling equation in step S32 is expressed as:

[0083] (26),

[0084] In formula (26),

[0085] ,

[0086] ,

[0087] ,

[0088] ,

[0089] In step S32, the forward equation (21) in step S31 is replaced by equation (26), and the target functional of step S32 is constructed using the same inversion method as step S31:

[0090] (27),

[0091] In formula (27), and Seismic data The variance of the noise in and the variance of the isotropic part of the model parameters, represents the regularization coefficient of the isotropic part of the model parameters, The low-frequency components of the isotropic part model parameters are represented. Finally, the iterative weighted least squares algorithm is used to solve the target functional in step S32 to obtain the model parameters M, μ, and ρ of the isotropic part.

[0092] The beneficial technical effects brought about by the present invention are:

[0093] The present invention utilizes the amplitude differences of seismic reflection data in different azimuths to invert the model parameters of the anisotropic part, and then uses the inverted model parameters of the anisotropic part as input to further invert the model parameters of the isotropic part. This method effectively reduces the pathological nature of multi-parameter simultaneous inversion and improves the stability and accuracy of the inversion. BRIEF DESCRIPTION OF THE DRAWINGS

[0094] Figure 1 This is a flow chart of an embodiment of the present invention.

[0095] Figure 2 Schematic diagram of the biaxial stress-induced anisotropic medium model.

[0096] Figure 3 Schematic diagram of the reflection coefficient changing with the observation azimuth angle.

[0097] Figure 4 Schematic diagram of the reflection coefficient changing with the incident angle.

[0098] Figure 5 is the model parameter Schematic diagram of the inversion results.

[0099] Figure 6 is the model parameter Schematic diagram of the inversion results.

[0100] Figure 7 is the model parameter Schematic diagram of the inversion results.

[0101] Figure 8 is the model parameter Schematic diagram of the inversion results.

[0102] Figure 9 is the model parameter Schematic diagram of the inversion results.

[0103] Figure 10 is the model parameter Schematic diagram of the inversion results.

[0104] Figure 11 is the model parameter Schematic diagram of the inversion results. DETAILED DESCRIPTION

[0105] The detailed description and technical contents of the present invention are described below with reference to the accompanying drawings. However, the drawings are only provided for reference and explanation and are not intended to limit the present invention.

[0106] like Figure 1 As shown, an embodiment of the present invention provides a method for azimuthal seismic inversion of biaxial stress-induced anisotropic media, comprising the following steps:

[0107] Step S1, parameterization of the stiffness tensor of the biaxial stress-induced anisotropic medium;

[0108] Step S2, deriving the PP wave reflection coefficient equation of the biaxial stress-induced anisotropic medium;

[0109] Step S3: azimuthal seismic inversion of biaxial stress-induced anisotropic media.

[0110] Step S1, the specific steps of parameterizing the stiffness tensor of the biaxial stress-induced anisotropic medium are as follows.

[0111] Considering that the unstressed rock is elastically isotropic and filled with randomly distributed and oriented microcracks, the elastic stiffness tensor of the unstressed rock can be expressed as follows according to the linear Hooke's law:

[0112] (1),

[0113] In formula (1), and represent the longitudinal wave modulus and shear modulus of unstressed rock, respectively.

[0114] like Figure 2 As shown, when isotropic rock is subjected to horizontal biaxial stress and When the stress acts, microcracks perpendicular to the stress direction will preferentially close, causing elastic anisotropy in the rock. Here, ψ represents the angle between the projection of the seismic wave propagation direction on the horizontal plane and the X1 axis, which is defined as the observation azimuth; θ represents the angle between the incident seismic wave propagation direction and the X3 axis, which is defined as the seismic wave incident angle.

[0115] According to Collet et al. (Geophysics, 2014, 79(5): C135-C145, doi: 10.1190 / geo2013-0311.1.), the equivalent flexibility tensor caused by microcrack closure due to horizontal biaxial stress can be expressed as:

[0116] (2),

[0117] In formula (2), Represents the equivalent flexibility matrix of rock without microcracks. and They represent the normal and tangential flexibility of microcracks in unstressed rock. The normalized stress can be expressed as and ,in, represents the typical closure pressure of the crack. and represents the compliance contribution of the initial microcracks, the tensor and Represents stress The flexibility contribution to crack closure, tensor and Represents stress The flexibility contribution that causes crack closure is expressed as follows:

[0118] (3),

[0119] (4),

[0120] (5),

[0121] In order to highlight the flexibility contribution of different parts, formula (2) is rearranged to obtain:

[0122] (6),

[0123] In formula (6), represents the flexibility tensor of the unstressed rock, and Represents horizontal stress 、 Compliance contribution to microcrack closure.

[0124] Furthermore, assuming that the biaxial horizontal stress is small and the anisotropy induced by the stress is weaker than that of unstressed rock, the equivalent stiffness tensor of the biaxial stress-induced anisotropic medium can be derived as:

[0125] (7),

[0126] Combining formulas (1)-(7), the stiffness tensor of biaxial stress-induced anisotropic media can be derived as:

[0127] (8),

[0128] From the symmetry of the stiffness tensor, it can be seen that the biaxial stress-induced anisotropic medium is approximately equivalent to the orthotropic medium. In order to simplify the stiffness tensor, the stress-induced anisotropy parameter is defined as follows:

[0129] (9),

[0130] (10),

[0131] Finally, the components of the elastic stiffness tensor in formula (8) can be expressed as:

[0132] (11),

[0133] In formula (11), , , .

[0134] Step S2, the specific steps of deriving the PP wave reflection coefficient equation of the biaxial stress-induced anisotropic medium are as follows.

[0135] According to Born scattering theory, the PP wave reflection coefficient equation can be expressed as:

[0136] (12),

[0137] In formula (12), represents the density of the background medium, is the angle of incidence, represents the perturbation of the stiffness tensor components, The specific expression is as follows:

[0138] (13),

[0139] In formula (13), represents the observation azimuth, Represents the longitudinal wave velocity of the background medium.

[0140] Consider the interface of two biaxial stress-induced anisotropic media with different properties. Combined with formula (11), the perturbation terms of each component of the stiffness tensor are derived. At the same time, assuming that the difference in properties of the media on both sides of the reflection interface is small, the effect of the , , , , , , , , , , , Proportional terms, we can get:

[0141] (14),

[0142] In formula (14), △ represents the difference in the properties of the media above and below the reflection interface.

[0143] Combining formulas (12)-(14), the PP wave reflection coefficient equation of biaxial stress-induced anisotropic media can be derived as follows:

[0144] (15),

[0145] In formula (15),

[0146] (16),

[0147] (17),

[0148] (18),

[0149] In formulas (16)-(18), , , ,

[0150] , ,

[0151] ,

[0152] .

[0153] Formula (15) is compared with the reflection coefficient equation of arbitrary weakly anisotropic media proposed by Pšenčík and Martins to verify the accuracy of the equation, as shown in Figure 3 and Figure 4 shown. Figure 3The curves of reflection coefficient changing with observation azimuth are shown. The dotted line represents the reflection coefficient curve calculated using the Pšenčík and Martins equations, and the solid line represents the reflection coefficient curve calculated using formula (15). Each curve corresponds to the incident angle θ = 0°, 10°, and 30° from top to bottom in the figure. Figure 4 The curve of reflection coefficient changing with incident angle is shown, where the dotted line represents the reflection coefficient curve calculated using Pšenčík and Martins equations, and the solid line represents the reflection coefficient curve calculated using formula (15). Each curve corresponds to the observation azimuth angle ψ = 90°, 45°, and 0° from top to bottom in the figure. Figure 3 and Figure 4 It can be seen that the PP wave reflection coefficient equation (15) proposed in the present invention for biaxial stress-induced anisotropic media has a high degree of agreement with the Pšenčík and Martins equations, proving the rationality of the derived equation.

[0154] Step S3, the specific steps of azimuthal seismic inversion of biaxial stress-induced anisotropic media are as follows.

[0155] Based on the azimuthal variation characteristics of the reflection coefficient of stress-induced anisotropic media, the model parameters of the anisotropic part are first inverted using the azimuthal amplitude difference method based on the Bayesian framework, considering the Cauchy prior information and low-frequency regularization constraints. Then, the inversion results are used as input to invert the model parameters of the isotropic part.

[0156] Step S31: According to formula (15), the azimuth difference reflection coefficient can be expressed as:

[0157] (19),

[0158] In formula (19), and represent two different observation azimuths, and

[0159] (20),

[0160] Considering the case of M reflection interfaces and L incident angles, based on the convolution model, the azimuth difference of seismic data can be expressed in matrix form as follows:

[0161] (twenty one),

[0162] In formula (21),

[0163] ,

[0164] The incident angle is And the observation azimuth is The earthquake trace;

[0165] ,

[0166] ,

[0167] ,

[0168] , , ,

[0169] , ,

[0170] Where W represents the wavelet matrix, Represents the Mth reflection interface. and They represent matrix transposition and generation of an M×M diagonal matrix respectively.

[0171] In the Bayesian theory framework, the posterior probability distribution of the model parameters is Proportional to the likelihood function of seismic data and the prior probability distribution of the model parameters ,Right now

[0172] (twenty two),

[0173] In formula (22), It is usually assumed to obey the Cauchy distribution, expressed as

[0174] (twenty three),

[0175] In formula (23), represents the variance of the model parameters.

[0176] Assuming that the likelihood function follows a Gaussian distribution, it can be expressed as:

[0177] (twenty four),

[0178] In formula (24), represents the variance of the noise in the seismic data. The background constraint of the model parameters can effectively improve the stability of the inversion. Therefore, combined with formulas (22)-(24), the target functional can be constructed as:

[0179] (25),

[0180] In formula (25), represents the regularization coefficient of the model parameters, Represents the low-frequency components of the model parameters.

[0181] The iterative weighted least squares algorithm is used to optimize and solve formula (25). After multiple iterations, the unknown model parameters can be obtained.

[0182] Step S32: The model parameters of the anisotropic part obtained by inversion in step S31 (assuming it is ) as input to invert the model parameters of the isotropic part.

[0183] The earthquake forward modeling equation in step S32 can be expressed as:

[0184] (26),

[0185] In formula (26),

[0186] ,

[0187] ,

[0188] ,

[0189] ,

[0190] In the inversion of step S32, the forward equation (21) in the inversion of step S31 is replaced by equation (26), and the inversion target functional of step S32 can be constructed using the same inversion method as step S31, and the following is obtained:

[0191] (27),

[0192] in, and Seismic data The variance of the noise in and the variance of the isotropic part of the model parameters, represents the regularization coefficient of the isotropic part of the model parameters, The low-frequency components of the isotropic model parameters are represented by the iterative weighted least squares algorithm. Finally, the target functional is solved to obtain the isotropic model parameters M, μ, and ρ.

[0193] Case study testing the inversion method:

[0194] An oil and gas reservoir in southwest China is used as an example to verify the inversion method proposed in the present invention. The oil and gas reservoir is a tight shale reservoir with low porosity and low permeability, and the structural fracture system is an important channel for oil and gas seepage. In the target work area, it is assumed that the oil-bearing reservoir can be regarded as an orthogonal medium induced by horizontal biaxial stress. Before being used in the present invention, the seismic traces collected on site will be preprocessed according to the data standard of azimuthal seismic inversion, including denoising, amplitude compensation, rotation angle domain, etc., and azimuthal seismic data will be generated. Then, azimuthal seismic inversion is carried out to estimate stress-induced anisotropy parameters and elastic parameters, and the results are as follows. Figure 5-Figure 11 As shown. It can be found that in the oil-bearing target reservoir (the location circled by the black oval frame), the inverted longitudinal modulus, shear modulus and density all show abnormally low values ​​( Figure 5-Figure 7 ). This indicates that fluid will reduce the elastic modulus of rock, which is consistent with the prediction of Gassmann's fluid substitution formula, and indirectly reflects the rationality of the inversion results of the isotropic model parameters. The inversion results of the stress-induced anisotropic parameters are obviously higher ( Figures 8-11 ), indicating that the target formation has significant stress-induced anisotropy. In addition, the correlation between the seven model parameters confirms the reliability of the inversion results.

[0195] The above embodiments are only preferred embodiments for fully illustrating the present invention, and the protection scope of the present invention is not limited thereto. Any equivalent substitution or modification made by those skilled in the art based on the present invention is within the protection scope of the present invention.

Claims

1. A method for azimuthal seismic inversion of biaxial stress-induced anisotropic media, characterized in that: It includes the following steps: Step S1, parameterization of the stiffness tensor of biaxial stress-induced anisotropic media; Step S2, deriving a PP wave reflection coefficient equation for a biaxial stress-induced anisotropic medium; Step S3: azimuthal seismic inversion of biaxial stress-induced anisotropic media.

2. A method for azimuthal seismic inversion of biaxial stress-induced anisotropic media according to claim 1, characterized in that: The step S1 comprises: Considering that the unstressed rock is elastically isotropic and full of randomly distributed and randomly oriented microcracks, the elastic stiffness tensor of the unstressed rock is expressed as: （1）， In formula (1), and represent the longitudinal wave modulus and shear modulus of unstressed rock respectively; The equivalent flexibility tensor caused by the closure of microcracks due to horizontal biaxial stress is expressed as: （2）， In formula (2), represents the equivalent flexibility matrix of the rock without microcracks, and They represent the normal flexibility and tangential flexibility of microcracks in unstressed rock, respectively. The normalized stress is expressed as and , The tensor represents the typical closure pressure of the fracture, and represents the compliance contribution of the initial microcrack, the tensor and Represents stress The flexibility contribution that causes crack closure, tensor and Represents stress The tensor expressions of the flexibility contribution that causes crack closure are: （3）， （4）， （5）， Rearranging formula (2) yields: （6）， In formula (6), represents the flexibility tensor of unstressed rock, and Represents horizontal stress , Compliance contribution that causes microcrack closure; Assuming biaxial stress-induced anisotropic media, the equivalent stiffness tensor is derived as: （7）， Combining formulas (1)-(7), the stiffness tensor of biaxial stress-induced anisotropic media is derived as follows: （8）， The biaxial stress-induced anisotropic medium is approximately equivalent to the orthotropic anisotropic medium, and the stress-induced anisotropic parameter is defined as: （9）， （10）， The components of the elastic stiffness tensor in formula (8) are expressed as: （11）， Formula (11) , , .

3. A method for azimuthal seismic inversion of biaxial stress-induced anisotropic media as claimed in claim 2, characterized in that: The step S2 comprises: The PP wave reflection coefficient equation is expressed as: （12）， In formula (12), represents the density of the background medium, is the angle of incidence, represents the perturbation of the stiffness tensor components, The expression is as follows: （13）， In formula (13), represents the observation azimuth, represents the longitudinal wave velocity of the background medium; Consider the interface of two biaxial stress-induced anisotropic media with different properties, and derive the perturbation terms of each component of the stiffness tensor by combining formula (11), ignoring the , , , , , , , , , , , Proportional terms, we get: （14）， In formula (14), △ represents the difference in the properties of the media above and below the reflection interface; Combining formulas (12)-(14), the PP wave reflection coefficient equation of biaxial stress-induced anisotropic media is derived as: （15）， In formula (15), （16）， （17）， （18）， Formula (16)-(18) , , , ， ， ， 。 4. A method for azimuthal seismic inversion of biaxial stress-induced anisotropic media as claimed in claim 3, characterized in that: The step S3 comprises: Step S31: According to formula (15), the azimuth difference reflection coefficient is expressed as: （19）， In formula (19), and Represent two different observation azimuths, and, （20）， Considering the case of M reflection interfaces and L incident angles, based on the convolution model, the azimuth difference of seismic data is expressed in matrix form as follows: （21）， In formula (21), ， The incident angle is And the observation azimuth is The seismic traces; ， ， ， ， ， ， ， ， Where W represents the wavelet matrix, represents the Mth reflection interface, and Respectively represent matrix transposition and generate an M×M diagonal matrix; In the Bayesian theory framework, the posterior probability distribution of the model parameters is Proportional to the likelihood function of the seismic data and the prior probability distribution of the model parameters : （22）， In formula (22), Assuming it obeys Cauchy distribution, it can be expressed as: （23）， In formula (23), represents the variance of the model parameters; Assuming that the likelihood function follows a Gaussian distribution, it can be expressed as: （24）， In formula (24), represents the variance of noise in seismic data; Combining formulas (22)-(24), the target functional of step S31 is constructed as: （25）， In formula (25), represents the regularization coefficient of the model parameters, Represents the low-frequency components of the model parameters; The iterative weighted least squares algorithm is used to optimize and solve formula (25). After multiple iterations, the unknown model parameters are obtained; Step S32: The model parameters of the anisotropic part obtained by inversion in step S31 are , as input, invert the model parameters of the isotropic part; The seismic forward modeling equation in step S32 is expressed as: （26）， In formula (26), ， ， ， ， In step S32, the forward equation (21) in step S31 is replaced by equation (26), and the target functional of step S32 is constructed using the same inversion method as step S31: （27）， In formula (27), and Seismic data The variance of the noise in and the variance of the isotropic part of the model parameters, represents the regularization coefficient of the isotropic part of the model parameters, The low-frequency components of the model parameters of the isotropic part are represented. Finally, the iterative weighted least squares algorithm is used to solve the target functional of step S32 to obtain the model parameters M, μ, and ρ of the isotropic part.

Images