Unsteady pre-stack stochastic seismic inversion methods, apparatus and equipment
By combining well logging and seismic data, a pre-stack attenuation forward model was constructed and stochastic inversion was performed, which solved the problems of numerical instability and poor lateral continuity in seismic inversion and achieved high-precision and high-resolution inversion results.
Patent Information
- Application Number
- CN202310808366.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-07-03
- Publication Date
- 2025-10-28
- Estimated Expiration
- 2043-07-03
AI Technical Summary
Existing seismic inversion methods suffer from numerical instability, poor lateral continuity, low accuracy and resolution in prestack inversion, especially inaccurate inversion results caused by amplitude attenuation and noise.
An initial inversion model is constructed based on well logging data to determine the prior information and probability distribution of the physical property parameters to be inverted. A pre-stack attenuation forward model is then constructed, and multiple results are generated through random inversion. Finally, the target inversion result is determined, and constraints are applied by combining the posterior information from well logging and seismic data.
The lateral continuity and resolution of the inversion results are improved, the accuracy of the inversion results is enhanced, and it meets the actual production needs.
Smart Images

Figure CN116819613B_ABST
Abstract
Description
Technical Field
[0001] This application relates to the field of petroleum exploration technology, and in particular to a non-steady-state pre-stack stochastic seismic inversion method, apparatus and equipment. Background Technology
[0002] This section is intended to provide background or context for the embodiments of the invention set forth in the claims. The description herein is not an admission that it is prior art simply because it is included in this section.
[0003] Seismic inversion methods play a crucial role in hydrocarbon prediction and reservoir characterization. They infer model parameters of the subsurface medium based on observed seismic data, established forward models, and prior information. Seismic inversion methods can include post-stack seismic inversion and pre-stack seismic inversion. Post-stack seismic inversion can only obtain a single physical property parameter, resulting in less than ideal performance in characterizing reservoir features and predicting fluids. Pre-stack seismic data contains more information on amplitude variations with offset, allowing pre-stack inversion to obtain multiple physical property parameters and better describe reservoir characteristics and predict fluids. The resolution of seismic inversion is affected by the accuracy of the seismic data. Due to the absorption and attenuation effect of the formation, seismic waves experience amplitude attenuation and phase distortion during propagation. This absorption and attenuation characteristic is characterized by the quality factor Q. Severe absorption and attenuation significantly impact the resolution of seismic data, leading to inaccurate inversion results. Therefore, existing methods compensate for seismic data before performing seismic inversion, such as inverse Q-filtering, which can effectively compensate for amplitude and correct phase, thereby eliminating the absorption and attenuation effect of the formation. However, due to the exponential decay of amplitude and the presence of random noise in the original seismic data, the compensated amplitude of attenuated seismic data exhibits instability. Furthermore, traditional compensation methods do not consider spatial continuity, making the compensation results susceptible to noise, resulting in poor lateral continuity, inaccurate inversion results, and low resolution.
[0004] Therefore, existing technologies cannot directly achieve high-resolution inversion of pre-stack physical property parameters of attenuated seismic data, and suffer from problems such as numerical instability, poor lateral continuity, and low accuracy and resolution. Summary of the Invention
[0005] The purpose of this application is to provide a non-steady-state pre-stack stochastic seismic inversion method, apparatus, and equipment to solve the problems of numerical instability, poor lateral continuity, and low accuracy and resolution in existing pre-stack inversion schemes.
[0006] To address the aforementioned technical problems, this specification provides a non-steady-state pre-stack stochastic seismic inversion method in its first aspect, comprising:
[0007] Based on the acquired well logging data, the initial inversion model is determined;
[0008] Based on the well logging data and the initial inversion model, the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data are determined.
[0009] Construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model;
[0010] Based on the probability distribution of the well logging data, the probability distribution of the seismic data, and the prior information of the physical property parameters to be inverted, the posterior information of the physical property parameters to be inverted is determined.
[0011] The posterior information is randomly inverted to generate multiple inversion results, and the target inversion result is determined based on the multiple inversion results.
[0012] In some embodiments, an initial inversion model is determined based on the acquired well logging data, including:
[0013] The well logging data is filtered and interpolated to obtain the initial inversion model.
[0014] In some embodiments, based on the well logging data and the initial inversion model, determining the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data includes:
[0015] Geostatistical analysis was performed on the well logging data to determine the geostatistical parameters of the physical property parameters to be inverted;
[0016] Based on the geostatistical parameters and the initial inversion model, the prior information of the physical property parameters to be inverted is determined;
[0017] Based on the prior information of the parameters to be inverted, and the relationship between the logging data and the physical property parameters to be inverted, the probability distribution of the logging data is determined.
[0018] In some embodiments, prior information of the physical property parameters to be inverted is determined based on the geostatistical parameters and the initial inversion model, including:
[0019] Based on the aforementioned geostatistical parameters, determine the spatial cross-correlation information of the physical property parameters to be inverted;
[0020] Based on the spatial cross-correlation information, the first covariance matrix of the physical property parameters to be inverted is determined;
[0021] The prior information is determined based on the initial inversion model and the first covariance matrix.
[0022] In some embodiments, a pre-stack attenuation forward model is constructed, and the probability distribution of attenuated seismic data is determined based on the pre-stack attenuation forward model, including:
[0023] The viscous filtering operator is integrated into the convolution model to construct a convolutional attenuation model as the pre-stack attenuation forward model.
[0024] Based on the convolutional attenuation model, the probability distribution of the seismic data is determined.
[0025] In some embodiments, the posterior information of the physical property parameters to be inverted includes the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and seismic data.
[0026] In some embodiments, the posterior information of the physical property parameter to be inverted is represented by the following formula:
[0027]
[0028] Where, μ m|(S,h) and Σ m|(S,h) Let μ represent the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and seismic data, respectively, in the posterior information. m Σ represents the mean in the prior information. m Let S represent the covariance matrix in the prior information, S represent the attenuated seismic data, h represent the well logging data, and C represent the covariance matrix in the prior information. (S,h) Σ represents the kernel matrix constrained by well logging and seismic data. e H represents the error between actual and simulated observation data, H represents the relationship operator between well logging data and the parameters to be inverted, and G represents the relationship operator between seismic data and the parameters to be inverted. T and G T Let H and G represent the operators obtained after transposing H and G respectively, and T represent the transpose symbol.
[0029] In some embodiments, random inversion is performed on the posterior information to generate multiple inversion results, and a target inversion result is determined based on the multiple inversion results, including:
[0030] The second covariance matrix is decomposed to obtain the target matrix;
[0031] The target matrix is randomly sampled using multiple random numbers to obtain multiple random sampling results;
[0032] The multiple inversion results are generated based on the multiple random sampling results and the mean in the posterior information;
[0033] The target inversion result is obtained by averaging the multiple inversion results.
[0034] The second aspect of this specification provides a non-steady-state pre-stack stochastic seismic inversion apparatus, comprising:
[0035] The initial model determination module is used to determine the initial inversion model based on the acquired well logging data;
[0036] The well logging data processing module is used to determine the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data based on the well logging data and the initial inversion model.
[0037] The seismic data processing module is used to construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model.
[0038] The posterior information generation module is used to determine the posterior information of the physical property parameter to be inverted based on the probability distribution of the well logging data, the probability distribution of the seismic data, and the prior information of the physical property parameter to be inverted.
[0039] The random inversion module is used to perform random inversion on the posterior information, generate multiple inversion results, and determine the target inversion result based on the multiple inversion results.
[0040] A third aspect of this specification provides a computer device, comprising: a memory and a processor, the processor and the memory being communicatively connected to each other, the memory storing computer instructions, and the processor executing the computer instructions to implement the steps of the method described in the first aspect.
[0041] A fourth aspect of this specification provides a computer-readable storage medium storing computer program instructions that, when executed, implement the steps of the method described in the first aspect.
[0042] The unsteady-state pre-stack stochastic seismic inversion method, apparatus, and equipment provided in this specification construct an initial inversion model using acquired well logging data. Based on the well logging data and the initial inversion model, prior information of the physical property parameters to be inverted and the probability distribution of the well logging data are determined. A pre-stack attenuation forward model is then constructed to determine the probability distribution of attenuated seismic data. Then, based on the determined probability distribution of attenuated seismic data, the probability distribution of well logging data, and the prior information, the posterior information of the physical property parameters to be inverted is determined. Subsequently, stochastic inversion can be performed based on the posterior information to obtain multiple inversion results, and the target inversion result of the pre-stack seismic inversion is determined based on these multiple results. This application uses attenuated seismic data for inversion, effectively avoiding the instability of amplitude compensation schemes used to compensate for attenuated seismic data during seismic inversion. Furthermore, this application uses both well logging data and seismic data as constraints to determine the posterior information of the physical property parameters to be inverted, enhancing the lateral continuity of the inversion results, improving the resolution and accuracy of the inversion results, and making the inversion results more consistent with actual production needs. Attached Figure Description
[0043] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this application. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.
[0044] Figure 1 The diagram shown is a flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0045] Figure 2 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0046] Figure 3 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0047] Figure 4 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0048] Figure 5 The diagram shown is a schematic representation of a true longitudinal wave velocity model provided in an embodiment of this application.
[0049] Figure 6 The diagram shown is a schematic representation of a true model of shear wave velocity provided in an embodiment of this application.
[0050] Figure 7 The diagram shown is a schematic representation of a density realism model provided in an embodiment of this application.
[0051] Figure 8 The diagram shown is a schematic diagram of an initial model of longitudinal wave velocity provided in an embodiment of this application;
[0052] Figure 9 The diagram shown is a schematic diagram of an initial model of shear wave velocity provided in an embodiment of this application;
[0053] Figure 10 The diagram shown is a schematic diagram of an initial density model provided in an embodiment of this application;
[0054] Figure 11 The diagram shown is a schematic representation of seismic data provided in an embodiment of this application.
[0055] Figure 12 The diagram shown is a schematic representation of attenuated seismic data provided in an embodiment of this application.
[0056] Figure 13The figure shown is a schematic diagram of the inversion result of a longitudinal wave velocity model provided in an embodiment of this application;
[0057] Figure 14 The figure shown is a schematic diagram of the inversion result of a shear wave velocity model provided in an embodiment of this application;
[0058] Figure 15 The figure shown is a schematic diagram of a density model inversion result provided in an embodiment of this application;
[0059] Figure 16 The image shown is a waveform diagram of a multiple random inversion result provided in an embodiment of this application;
[0060] Figure 17 The image shown is a vibrational spectrum provided in an embodiment of this application;
[0061] Figure 18 The diagram shown is a structural schematic of a non-steady-state pre-stack stochastic seismic inversion device provided in an embodiment of this application.
[0062] Figure 19 The diagram shown is a schematic block diagram of an electronic device provided in an embodiment of this application. Detailed Implementation
[0063] To enable those skilled in the art to better understand the technical solutions in this application, the technical solutions in the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of this application, and not all of the embodiments. Based on the embodiments in this application, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of this application.
[0064] As mentioned earlier, existing methods compensate for seismic data during seismic inversion to eliminate the absorption and attenuation effects of the formation. However, since amplitude decays exponentially and the original seismic data usually contains random noise, the compensated amplitude of the attenuated seismic data is unstable. Furthermore, traditional compensation methods do not consider spatial continuity, making the compensation results susceptible to noise, resulting in poor lateral continuity, inaccurate inversion results, and low resolution.
[0065] Considering the impact of amplitude compensation on the accuracy and resolution of seismic inversion, this application proposes an unsteady pre-stack stochastic seismic inversion method, apparatus, and equipment to address the aforementioned issues. The method includes: determining an initial inversion model based on acquired well logging data; determining prior information of the physical property parameters to be inverted and the probability distribution of the well logging data based on the well logging data and the initial inversion model; constructing a pre-stack attenuation forward model and determining the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model; determining posterior information of the physical property parameters to be inverted based on the probability distribution of the well logging data, the probability distribution of the seismic data, and the prior information of the physical property parameters to be inverted; performing stochastic inversion on the posterior information to generate multiple inversion results, and determining a target inversion result based on the multiple inversion results.
[0066] This application utilizes attenuated seismic data for inversion, integrating the traditional two-step inversion method (compensation followed by inversion) into a single-step inversion method that requires no compensation. This effectively avoids the instability of amplitude compensation schemes used to compensate for attenuated seismic data during seismic inversion. Furthermore, this application uses well logging data and seismic data as constraints, and within a Bayesian theoretical framework, determines the posterior information of the physical property parameters to be inverted. Combined with stochastic inversion methods, this enhances the lateral continuity of the inversion results while improving their resolution and accuracy, making the inversion results more consistent with actual production needs.
[0067] The method provided in this application can be executed by an electronic device, which is an electronic device with data computing, processing, and storage capabilities. This electronic device can be a terminal such as a personal computer (PC), tablet computer, smartphone, wearable device, or intelligent robot; or it can be a server. The server can be an independent physical server, a server cluster or distributed system composed of multiple physical servers, or a cloud server providing cloud computing services.
[0068] The following section, in conjunction with the accompanying drawings, introduces the non-steady-state pre-stack stochastic seismic inversion method provided in the embodiments of this application.
[0069] Figure 1 The diagram shown is a flowchart illustrating a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application. Figure 1 As shown, the method includes:
[0070] S101: Determine the initial inversion model based on the acquired well logging data.
[0071] Well logging data is understood to be data that characterizes geological structures and formation changes through observation. It can include both high-frequency and low-frequency information, and noise interference can be reduced by filtering the data.
[0072] In some embodiments, step S101, based on the acquired well logging data, determines the initial inversion model, including:
[0073] The well logging data is filtered and interpolated to obtain the initial inversion model.
[0074] It is understandable that, since the noise signal is a high-frequency signal, low-pass filtering is performed on the logging data to obtain denoised logging data. Furthermore, interpolation and extrapolation are performed on the low-frequency logging data to obtain the low-frequency model, which is the initial model required for inversion (i.e., the initial inversion model in step S101).
[0075] It is understandable that the initial inversion model can be used as the mean of the physical property parameters to be inverted (which can be expressed as μ). m The prior information of the physical property parameters to be inverted can be further solved through this initial inversion model.
[0076] In some embodiments, the physical property parameters to be inverted may include one or more parameters such as P-wave velocity, S-wave velocity, density, and porosity. Furthermore, the initial inversion model may include one or more of the following: an initial P-wave velocity model, an initial S-wave velocity model, an initial density model, and an initial porosity model. It is understood that the above examples of the physical property parameters to be inverted and the initial inversion model are merely one embodiment of this application, and other physical property parameters describing reservoir characteristics and fluid predictions are also within the scope of protection of this application. The above examples are not intended to limit the scope of protection of this application.
[0077] S102: Based on the well logging data and the initial inversion model, determine the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data.
[0078] It is understandable that the prior information can be a prior probability distribution, which can be characterized by the mean and covariance matrix of the physical property parameters to be inverted.
[0079] Figure 2 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0080] like Figure 2 As shown, in some embodiments, based on the well logging data and the initial inversion model, the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data are determined, including:
[0081] S201: Perform geostatistical analysis on the well logging data to determine the geostatistical parameters of the physical property parameters to be inverted.
[0082] It is understandable that geostatistical analysis can simulate the physical property parameters to be inverted and obtain the corresponding geostatistical parameters.
[0083] In some embodiments, the physical property parameters to be inverted are continuous parameters, and a two-point geostatistical algorithm can be used to analyze the well logging data to simulate continuous physical property parameters such as P-wave velocity, S-wave velocity, and density. In other embodiments, the physical property parameters to be inverted are discrete parameters, and a multi-point geostatistical algorithm can be used to analyze the well logging data to simulate discrete physical property parameters such as lithofacies.
[0084] In some embodiments, geostatistical parameters may include range, sill value, etc., used to determine the spatial cross-correlation information of the physical property parameters to be inverted.
[0085] S202: Based on the geostatistical parameters and the initial inversion model, determine the prior information of the physical property parameters to be inverted.
[0086] In some embodiments, when determining prior information, the first covariance matrix (which can be represented as Σ) of the prior information of the physical property parameters to be inverted can be determined based on geostatistical parameters. m Furthermore, since the initial inversion model can be used as the mean value of the physical property parameters to be inverted is μ... m Then, based on the first covariance matrix Σ m The sum and mean are μ m This allows us to determine the prior information of the physical property parameters to be inverted.
[0087] Figure 3 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0088] like Figure 3 As shown, in some embodiments, in step S202, based on the geostatistical parameters and the initial inversion model, the prior information of the physical property parameters to be inverted is determined, including:
[0089] S301: Based on the geostatistical parameters, determine the spatial cross-correlation information of the physical property parameters to be inverted.
[0090] It can be understood that spatial cross-correlation information is the similarity of the spatial location values of the physical property parameter to be inverted. In some embodiments, the spatial cross-correlation information of the physical property parameter to be inverted can be characterized by a variogram function. Therefore, step S301 above may specifically include:
[0091] Based on geostatistical parameters—range and sill values—the variogram function of the physical property parameters to be inverted is determined. The variogram function can be determined using the following formula:
[0092]
[0093] Where h represents the spatial distance between the physical property parameters to be inverted, a represents the range, C represents the sill value, and γ(h) represents the variogram function.
[0094] Based on the above formula (3) and the determined geostatistical parameters, the variogram function between the physical property parameters to be inverted can be determined to characterize the spatial cross-correlation information of the physical property parameters to be inverted. The variogram function can better describe the distribution of the physical property parameters to be inverted.
[0095] S302: Based on the spatial cross-correlation information, determine the first covariance matrix of the physical property parameters to be inverted.
[0096] It is understood that there is a certain relationship between the spatial cross-correlation information and the covariance matrix of the physical property parameters to be inverted. Based on this relationship and the spatial cross-correlation information determined in formula S301, the first covariance matrix of the physical property parameters to be inverted can be calculated. The first covariance matrix can be understood as the covariance matrix of the prior information of the physical property parameters to be inverted.
[0097] In some embodiments, the first covariance matrix can be determined by the following formula:
[0098] C(h)=C-γ(h) Formula (2)
[0099] Where C represents the sill value of the geological statistical parameter, γ(h) represents the variogram function, i.e., the spatial cross-correlation information determined in step S301, and C(h) represents the covariance. Based on the covariance obtained from formula (2), the first covariance matrix Σ of the physical property parameters to be inverted can be determined. m .
[0100] S303: Determine the prior information based on the initial inversion model and the first covariance matrix.
[0101] Specifically, in some embodiments, the initial inversion model is the mean μ of the physical property parameters to be inverted. m Assume that the physical property parameters to be inverted follow a mean of μ. m The first covariance matrix is Σ m The multivariate Gaussian distribution is based on the first covariance matrix being Σ. m The sum and mean are μ m This allows us to determine the prior distribution function of the physical property parameters to be inverted. The prior distribution function can be characterized as: m:N(μ m ,Σm ). Here, m can represent the physical property parameter to be inverted.
[0102] S203: Based on the prior information of the parameters to be inverted, and the relationship between the logging data and the physical property parameters to be inverted, determine the probability distribution of the logging data.
[0103] It is understandable that there is a certain linear relationship between well logging data and the physical property parameters to be inverted. Through this linear relationship, the probability distribution of well logging data can be determined, thereby constraining the posterior information by well logging data and attenuated seismic data.
[0104] In some embodiments, the relationship between well logging data and the physical property parameters to be inverted can be represented by the formula h = Hm, where h represents well logging data, m represents the physical property parameters to be inverted, and H represents the forward modeling operator.
[0105] Furthermore, since the physical property parameters to be inverted follow a mean of μ m The first covariance matrix is Σ m The data follows a multivariate Gaussian distribution. Combined with the relationship between the logging data and the physical property parameters to be inverted, the logging data follows a mean of Hμ. m The covariance is HΣ m H T The probability distribution of well logging data can be represented by a multivariate Gaussian distribution as h:N(Hμ m ,HΣ m H T ).
[0106] S103: Construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model.
[0107] It is understandable that the pre-stack attenuation forward model is a pre-stack forward model constructed to characterize the stratigraphic filtering effect, taking into account the impact of amplitude compensation on the accuracy and resolution of seismic inversion. The pre-stack attenuation forward model can characterize the relationship between seismic data and seismic wavelets and reflection coefficients.
[0108] It is understandable that attenuated seismic data is seismic data acquired after being filtered by the strata.
[0109] In some embodiments, step S103 involves constructing a pre-stack attenuation forward model and determining the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model, including:
[0110] A viscous filtering operator is integrated into the convolution model to construct a convolutional attenuation model as the pre-stack attenuation forward model; based on the convolutional attenuation model, the probability distribution of the seismic data is determined.
[0111] In some embodiments, the convolution model can be represented by the following formula:
[0112]
[0113] Where, ω k The angular frequency is represented by i, the imaginary unit is represented by j, the time sampling point is represented by j, and the frequency sampling point is represented by k. and Representing the earthquake data s(τ) j ) and seismic wavelet w(τ) j In the frequency domain, r represents the reflection coefficient, and τ j This indicates the travel time of seismic waves.
[0114] Furthermore, in some embodiments, the convolutional decay model can be expressed by the following formula:
[0115]
[0116] Among them, Q ej Representing time 0-τ j The equivalent Q value between them Q j′ Representing time τ j′ The Q value at point ω, where Δt represents the time sampling interval, and α(ω) k ,τ j Q ej ) represents the viscous filtering operator, which represents the filtering effect of the formation.
[0117] Where α(ω) k ,τ j Q ej This can be expressed by the following formula:
[0118]
[0119] In formula (5), ω q The Nyquist angular frequency, a j Indicates the attenuation factor.
[0120] In some embodiments, the convolutional attenuation model in formula (4) can be written in matrix-vector form to obtain the following formula:
[0121]
[0122] Among them, matrix Represents the seismic gathers in the frequency domain, with matrix R representing the reflection coefficients. The matrix represents the filtering effect of formation Q. This represents the frequency domain seismic wavelet.
[0123] Furthermore, an inverse Fourier transform is performed on the above formula (6), that is, both sides of formula (6) are multiplied by the inverse Fourier transform factor e. iωt This allows us to obtain the time-domain pre-stack forward modeling equations for attenuated seismic data. Specifically, it can be expressed by the following formula:
[0124] S=WΛADm formula (7)
[0125] Let G = WΛAd, where G represents the forward operator, and we can obtain the following formula:
[0126] S=Gm Formula (8)
[0127] Formula (8) can be used as the pre-stack attenuation forward model constructed in step S103.
[0128] Assume that the attenuated seismic data follows a multivariate Gaussian distribution S:N(Gμ) m ,GΣ m G T +Σ e The mean value of the attenuated seismic data is Gμ. m The covariance is GΣ m G T +Σ e .
[0129] S104: Based on the probability distribution of the well logging data, the probability distribution of the seismic data, and the prior information of the physical property parameters to be inverted, determine the posterior information of the physical property parameters to be inverted.
[0130] Based on linear Bayesian theory, the posterior probability of the physical property parameter to be inverted can be determined based on the likelihood function and the prior probability. Specifically, the posterior probability can be expressed by the following formula:
[0131] p(α|β)∝p(β|α)p(α) Formula (9)
[0132] Where p(α|β) represents the posterior probability of the physical property parameter to be inverted, p(β|α) represents the likelihood function, p(α) represents the prior probability of the physical property parameter to be inverted, and ∝ indicates that the posterior probability is proportional to the product of the likelihood function and the prior probability.
[0133] In some embodiments, the posterior information of the physical property parameters to be inverted includes the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and seismic data.
[0134] It is understandable that, in order to improve the lateral continuity of the inversion results, the posterior probability in formula (9) can be constrained by combining well logging data and attenuated seismic data based on linear Bayesian theory.
[0135] Specifically, in some embodiments, the prior probability of the physical property parameter m to be inverted follows a mean of μ. m The covariance is Σ m The multivariate Gaussian distribution can be represented as m:N(μ m ,Σ m Well logging data follows a mean of Hμ m The covariance is HΣ m H T The multivariate Gaussian distribution can be represented as h:N(Hμ m ,HΣ m H T The attenuated seismic data follows a mean of Gμ. m The covariance is GΣ m G T +Σ e The multivariate Gaussian distribution can be represented as S:N(Gμ m ,GΣ m G T +Σ e ).
[0136] The posterior information of the physical property parameters to be inverted is expressed by the following formula:
[0137]
[0138] Where, μ m|(S,h) and Σ m|(S,h) Let μ represent the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and seismic data, respectively, in the posterior information. m Σ represents the mean in the prior information. m Let S represent the covariance matrix in the prior information, S represent the attenuated seismic data, h represent the well logging data, and C represent the covariance matrix in the prior information. (S,h) Σ represents the kernel matrix constrained by well logging and seismic data. e H represents the error between actual and simulated observation data, H represents the relationship operator between well logging data and the parameters to be inverted, and G represents the relationship operator between seismic data and the parameters to be inverted. T and G T Let H and G represent the operators obtained after transposing H and G respectively, and T represent the transpose symbol.
[0139] It is understandable that the mean μ m|(S,h) The maximum value of the posterior probability distribution of the physical property parameters to be inverted represents the optimal solution of the inversion result, which is a deterministic result. The second covariance matrix Σ m|(S,h)This can be used for subsequent stochastic simulations to obtain multiple reasonable inversion results. Furthermore, by combining multiple inversion results obtained from the stochastic inversion of the second covariance matrix with the mean of the optimal solution of the inversion results, multiple reasonable inversion results can be obtained. Based on these multiple reasonable inversion results, the target inversion result of the pre-stack seismic inversion can be determined.
[0140] S105: Perform random inversion on the posterior information to generate multiple inversion results, and determine the target inversion result based on the multiple inversion results.
[0141] Specifically, the stochastic inversion of posterior information can be the stochastic inversion of the second covariance matrix in the posterior information.
[0142] Figure 4 The diagram shown is a partial flowchart of a non-steady-state pre-stack stochastic seismic inversion method provided in an embodiment of this application.
[0143] like Figure 4 As shown, in some embodiments, random inversion is performed on the posterior information to generate multiple inversion results, and a target inversion result is determined based on the multiple inversion results, including:
[0144] S401: Decompose the second covariance matrix to obtain the target matrix.
[0145] It is understandable that performing a random inversion on the second covariance matrix directly would be computationally intensive. Therefore, in order to reduce the computational burden of the inversion, the second covariance matrix can be decomposed to simplify the calculation.
[0146] Specifically, the second covariance matrix can be decomposed into the product of an upper triangular matrix and a lower triangular matrix using the LU decomposition technique. That is, the second covariance matrix can be expressed by the following formula:
[0147] Σ m|(S,h) =LU Formula (12)
[0148] Here, matrix L represents an upper triangular matrix, matrix U represents a lower triangular matrix, and L = U T .
[0149] It is understood that the target matrix in step S401 can be the upper triangular matrix in formula (12) or the lower triangular matrix in formula (12), and this application does not limit it.
[0150] S402: Randomly sample the target matrix using multiple random numbers to obtain multiple random sampling results.
[0151] Specifically, the target matrix can be randomly sampled using multiple random numbers r. Taking an upper triangular matrix as an example, the sampling result in step S402 can be represented as Lr, where L represents the upper triangular matrix, and r represents the random numbers used for random sampling, which are random numbers between 0 and 1 that follow a uniform distribution.
[0152] S403: Based on the multiple random sampling results and the mean in the posterior information, generate the multiple inversion results. Specifically, the inversion results can be expressed by the following formula:
[0153] I m =μ m|(S,h) +Lr formula (13)
[0154] Among them, I m Indicates the inversion result, μ m|(S,h) Let L represent the mean in the posterior information, L represent the upper triangular matrix (i.e., the target matrix), and r represent the random number used for random sampling.
[0155] S404: Take the average value of the multiple inversion results to obtain the target inversion result.
[0156] The effects of the unsteady pre-stack stochastic seismic inversion method in this application will be described below with reference to the accompanying drawings and specific embodiments.
[0157] In this embodiment, the real model includes a true P-wave velocity model, a true S-wave velocity model, and a true density model; the initial inversion model includes an initial P-wave velocity model, an initial S-wave velocity model, and an initial density model.
[0158] Figures 5 to 7 The figures shown are schematic diagrams of the true P-wave velocity model, the true S-wave velocity model, and the true density model provided in the embodiments of this application.
[0159] Depend on Figures 5-7 As can be seen, the real model includes 737 seismic traces with a longitudinal time of 750ms, and serves as the real model for reference.
[0160] Figures 8-10 The figures shown are schematic diagrams of the initial models for longitudinal wave velocity, transverse wave velocity, and density provided in the embodiments of this application.
[0161] Figure 11 The diagram shown is a schematic representation of seismic data provided in an embodiment of this application. Specifically, Figure 11 The earthquake data in the data are pre-stack earthquake data without attenuation obtained based on formula (3), and can be used as reference data.
[0162] Figure 12The diagram shown is a schematic representation of attenuated seismic data provided in an embodiment of this application. Specifically, Figure 12 The earthquake data shown are simulated observation earthquake data obtained based on formula (8).
[0163] Figures 13 to 15 This is a schematic diagram of the P-wave velocity model, S-wave velocity model, and density model obtained by the non-steady-state pre-stack stochastic seismic inversion method provided in the embodiments of this application.
[0164] Figure 16 The image shown is a waveform diagram of a multiple random inversion result provided in an embodiment of this application. Specifically, Figure 16 From Figures 13-15 The provided target inversion results include multiple random simulation results of the single-channel inversion results at the well location (CDP=100), where the solid black line represents the true model. Figures 5-7 In the model, the gray line represents the results of multiple random inversions, and the dashed black line represents the average of the results of multiple random simulations, which is used as the final inversion result. Figure 16 The three figures in the image, from left to right, represent the target inversion results waveforms for longitudinal wave velocity, transverse wave velocity, and density, respectively.
[0165] Figure 17 The image shown is a vibrational spectrum provided in an embodiment of this application. Wherein, Figure 17 The corresponding non-attenuated seismic data in Figure 11 The single-channel seismic data with channel number 100 has the following attenuation data: Figure 12 The middle channel number corresponds to the single-channel seismic data of 100, and the composite record of the inversion result is as follows: Figure 17 The mean values of the inversion results (P-wave velocity Vp, S-wave velocity Vs, density ρ) represented by the dashed and solid lines are used to obtain the synthetic seismic record using formula (3).
[0166] Depend on Figures 5-7 as well as Figures 13-15 It can be seen that the inversion results of the physical property parameters obtained in this application have a high degree of agreement with the actual model. Figure 17 It can be seen that the unattenuated seismic records and the composite records obtained from the inversion results have a high degree of agreement.
[0167] Therefore, the non-steady-state pre-stack stochastic seismic inversion method provided in this application can achieve good amplitude compensation while obtaining high-resolution and high-precision inversion results.
[0168] This application also provides an unsteady pre-stack stochastic seismic inversion device. Figure 18 The diagram shown is a structural schematic of the unsteady pre-stack stochastic seismic inversion device 1800 provided in an embodiment of this application. Figure 18As shown, the device includes:
[0169] The initial model determination module 1801 is used to determine the initial inversion model based on the acquired well logging data;
[0170] The well logging data processing module 1802 is used to determine the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data based on the well logging data and the initial inversion model.
[0171] The seismic data processing module 1803 is used to construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model.
[0172] The posterior information generation module 1804 is used to determine the posterior information of the physical property parameter to be inverted based on the probability distribution of the well logging data, the probability distribution of the seismic data, and the prior information of the physical property parameter to be inverted.
[0173] The random inversion module 1805 is used to perform random inversion on the posterior information, generate multiple inversion results, and determine the target inversion result based on the multiple inversion results.
[0174] The descriptions and functions of the above modules can be understood by referring to the section on unsteady pre-stack stochastic seismic inversion methods, and will not be repeated here.
[0175] This application also provides an electronic device, such as... Figure 19 As shown, the electronic device may include a processor 1901 and a memory 1902, wherein the processor 1901 and the memory 1902 may be connected via a bus or other means. Figure 19 Taking the example of a connection between China and Israel via a bus.
[0176] Processor 1901 can be a Central Processing Unit (CPU). Processor 1901 can also be other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, or combinations of the above types of chips.
[0177] Memory 1902, as a non-transitory computer-readable storage medium, can be used to store non-transitory software programs, non-transitory computer-executable programs, and modules, such as the program instructions / modules corresponding to the non-steady-state pre-stack stochastic seismic inversion method in this embodiment of the invention (e.g., Figure 18 The module shown includes an initial model determination module 1801, a well logging data processing module 1802, a seismic data processing module 1803, a posterior information generation module 1804, and a stochastic inversion module 1805. The processor 1901 executes various functional applications and data processing by running non-transient software programs, instructions, and modules stored in the memory 1902, thereby implementing the non-steady-state pre-stack stochastic seismic inversion method described in the above method embodiment.
[0178] The memory 1902 may include a program storage area and a data storage area. The program storage area may store the operating system and applications required for at least one function; the data storage area may store data created by the processor 1901, etc. Furthermore, the memory 1902 may include high-speed random access memory and may also include non-transitory memory, such as at least one disk storage device, flash memory device, or other non-transitory solid-state storage device. In some embodiments, the memory 1902 may optionally include memory remotely located relative to the processor 1901, and these remote memories may be connected to the processor 1901 via a network. Examples of such networks include, but are not limited to, the Internet, corporate intranets, local area networks, mobile communication networks, and combinations thereof.
[0179] The one or more modules are stored in the memory 1902, and when executed by the processor 1901, they perform the following: Figure 1 The non-steady-state pre-stack stochastic seismic inversion method in the illustrated embodiment.
[0180] The specific details of the aforementioned electronic device can be understood by referring to the relevant descriptions and effects in the above method embodiments, and will not be repeated here.
[0181] This specification also provides a computer storage medium storing computer program instructions that, when executed, implement the steps of the above-described unsteady pre-stack stochastic seismic inversion method.
[0182] Those skilled in the art will understand that all or part of the processes in the methods of the above embodiments can be implemented by a computer program instructing related hardware. The program can be stored in a computer-readable storage medium, and when executed, it can include the processes of the embodiments of the above methods. The storage medium can be a magnetic disk, optical disk, read-only memory (ROM), random access memory (RAM), flash memory, hard disk drive (HDD), or solid-state drive (SSD), etc.; the storage medium can also include combinations of the above types of memory.
[0183] The various embodiments in this specification are described in a progressive manner. For the same or similar parts between the various embodiments, please refer to each other. The focus of each embodiment is to describe the differences from other embodiments.
[0184] The systems, devices, modules, or units described in the above embodiments can be implemented by computer chips or entities, or by products with certain functions.
[0185] For ease of description, the above devices are described separately by function as various units. Of course, in implementing this application, the functions of each unit can be implemented in one or more software and / or hardware.
[0186] As can be seen from the above description of the embodiments, those skilled in the art can clearly understand that this application can be implemented by means of software plus necessary general-purpose hardware platforms. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product can be stored in a storage medium, such as ROM / RAM, magnetic disk, optical disk, etc., and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute certain parts of the methods of various embodiments of this application.
[0187] This application can be used in a wide variety of general-purpose or special-purpose computer system environments or configurations. For example: personal computers, server computers, handheld or portable devices, tablet devices, multiprocessor systems, microprocessor-based systems, set-top boxes, programmable consumer electronics devices, network PCs, minicomputers, mainframe computers, distributed computing environments including any of the above systems or devices, etc.
[0188] This application can be described in the general context of computer-executable instructions, such as program modules, that are executed by a computer. Generally, program modules include routines, programs, objects, components, data structures, etc., that perform a specific task or implement a specific abstract data type. This application can also be practiced in distributed computing environments where tasks are performed by remote processing devices connected via a communication network. In distributed computing environments, program modules can reside in local and remote computer storage media, including storage devices.
[0189] Although this application has been described through embodiments, those skilled in the art will know that this application has many modifications and variations without departing from the spirit of this application, and it is intended that the appended claims cover such modifications and variations without departing from the spirit of this application.
Claims
1. A non-steady-state pre-stack stochastic seismic inversion method, characterized in that, include: Based on the acquired well logging data, the initial inversion model is determined; Based on the well logging data and the initial inversion model, the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data are determined. Construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model; Based on the probability distribution of the well logging data, the probability distribution of the attenuated seismic data, and the prior information of the physical property parameters to be inverted, the posterior information of the physical property parameters to be inverted is determined. The posterior information of the physical property parameters to be inverted includes the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and the attenuated seismic data. The posterior probability of the physical property parameters to be inverted under the constraints of the well logging data and the attenuated seismic data follows a probability with a mean of μ. m|(S,h) The second covariance matrix is Σ m|(S,h) The multivariate Gaussian distribution; The posterior information is randomly inverted to generate multiple inversion results, and the target inversion result is determined based on the multiple inversion results; Constructing a pre-stack attenuation forward model and determining the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model, including: The viscous filtering operator is integrated into the convolution model to construct a convolutional attenuation model as the pre-stack attenuation forward model. Based on the convolutional attenuation model, the probability distribution of attenuated seismic data is determined; the attenuated seismic data follows a multivariate Gaussian distribution S:N(Gμ m ,GΣ m G T +Σ e The mean value of the attenuated seismic data is Gμ. m The covariance is GΣ m G T +Σ e In this context, the forward modeling operator G = WΛAd, where W represents the seismic wavelet and Λ represents the filtering effect of the formation Q; μ m Σ represents the mean in the prior information. m Let S represent the covariance matrix in the prior information, and let S represent the attenuated seismic data. e This indicates the error between actual observation data and simulated observation data.
2. The non-steady-state pre-stack stochastic seismic inversion method according to claim 1, characterized in that, Based on the acquired well logging data, an initial inversion model is determined, including: The well logging data is filtered and interpolated to obtain the initial inversion model.
3. The non-steady-state pre-stack stochastic seismic inversion method according to claim 1, characterized in that, Based on the well logging data and the initial inversion model, the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data are determined, including: Geostatistical analysis was performed on the well logging data to determine the geostatistical parameters of the physical property parameters to be inverted; Based on the geostatistical parameters and the initial inversion model, the prior information of the physical property parameters to be inverted is determined; Based on the prior information of the physical property parameters to be inverted, and the relationship between the well logging data and the physical property parameters to be inverted, the probability distribution of the well logging data is determined.
4. The unsteady pre-stack stochastic seismic inversion method according to claim 3, characterized in that, Based on the geostatistical parameters and the initial inversion model, the prior information of the physical property parameters to be inverted is determined, including: Based on the aforementioned geostatistical parameters, determine the spatial cross-correlation information of the physical property parameters to be inverted; Based on the spatial cross-correlation information, the first covariance matrix of the physical property parameters to be inverted is determined; The prior information is determined based on the initial inversion model and the first covariance matrix.
5. The unsteady pre-stack stochastic seismic inversion method according to claim 1, characterized in that, The posterior information of the physical property parameters to be inverted is expressed by the following formula: Where, μ m|(S,h) and Σ m|(S,h) Let μ represent the mean and second covariance matrix of the physical property parameters to be inverted under the constraints of the well logging data and seismic data, respectively, in the posterior information. m Σ represents the mean in the prior information. m Let S represent the covariance matrix in the prior information, S represent the attenuated seismic data, h represent the well logging data, and C represent the covariance matrix in the prior information. (S,h) Σ represents the kernel matrix constrained by well logging and seismic data. e H represents the error between actual and simulated observation data, H represents the relationship operator between well logging data and the parameters to be inverted, and G represents the relationship operator between seismic data and the parameters to be inverted. T and G T Let H and G represent the operators obtained after transposing H and G respectively, and T represent the transpose symbol.
6. The unsteady pre-stack stochastic seismic inversion method according to claim 1, characterized in that, Random inversion is performed on the posterior information to generate multiple inversion results, and a target inversion result is determined based on the multiple inversion results, including: The second covariance matrix is decomposed to obtain the target matrix; The target matrix is randomly sampled using multiple random numbers to obtain multiple random sampling results; The multiple inversion results are generated based on the multiple random sampling results and the mean in the posterior information; The target inversion result is obtained by averaging the multiple inversion results.
7. A non-steady-state pre-stack stochastic seismic inversion device, characterized in that, include: The initial model determination module is used to determine the initial inversion model based on the acquired well logging data; The well logging data processing module is used to determine the prior information of the physical property parameters to be inverted and the probability distribution of the well logging data based on the well logging data and the initial inversion model. The seismic data processing module is used to construct a pre-stack attenuation forward model and determine the probability distribution of attenuated seismic data based on the pre-stack attenuation forward model. The posterior information generation module is used to determine the posterior information of the physical property parameter to be inverted based on the probability distribution of the well logging data, the probability distribution of the attenuated seismic data, and the prior information of the physical property parameter to be inverted. The posterior information of the physical property parameter to be inverted includes the mean and second covariance matrix of the physical property parameter to be inverted under the constraints of the well logging data and the attenuated seismic data. The posterior probability of the physical property parameter to be inverted under the constraints of the well logging data and the attenuated seismic data follows a probability with a mean of μ. m|(S,h) The second covariance matrix is Σ m|(S,h) The multivariate Gaussian distribution; The random inversion module is used to perform random inversion on the posterior information, generate multiple inversion results, and determine the target inversion result based on the multiple inversion results; The seismic data processing module is specifically used for: The viscous filtering operator is integrated into the convolution model to construct a convolutional attenuation model as the pre-stack attenuation forward model. Based on the convolutional attenuation model, the probability distribution of attenuated seismic data is determined; the attenuated seismic data follows a multivariate Gaussian distribution S:N(Gμ m ,GΣ m G T +Σ e The mean value of the attenuated seismic data is Gμ. m The covariance is GΣ m G T +Σ e In this context, the forward modeling operator G = WΛAd, where W represents the seismic wavelet and Λ represents the filtering effect of the formation Q; μ m Σ represents the mean in the prior information. m Let S represent the covariance matrix in the prior information, and let S represent the well logging data. e This indicates the error between actual observation data and simulated observation data.
8. A computer device, characterized in that, include: A memory and a processor, the processor and the memory being communicatively connected to each other, the memory storing computer instructions, the processor executing the computer instructions to implement the steps of the method according to any one of claims 1 to 6.
9. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores computer program instructions that, when executed, implement the steps of the method according to any one of claims 1 to 6.
Citation Information
Patent Citations
Seismic random inversion method and device based on multi-point geostatistical prior information
CN110031896A
Viscoelastic medium seismic inversion method based on Zoeppritz equation
CN113050162A
Horizon-constrained multivariate Gaussian fast prestack random inversion method
CN115951410A