Unsupervised brittle parameter inversion method, device and medium

Abstract

The present invention belongs to the field of seismic inversion in oil and gas exploration, and specifically discloses an unsupervised brittle parameter inversion method, device and medium. This method is mainly applied to predict brittle parameters and optimize reservoir sections with good fracturability. It integrates a DL framework and an inversion method based on brittle parameter inversion equations, providing a solution strategy without training labels for DL pre-stack seismic inversion work. At the same time, it improves the inversion accuracy and efficiency and is applicable to large-scale data processing. The present invention uses a theoretical data and an actual data of a work area for experimental analysis, further verifying the inversion accuracy and efficiency of this method and verifying the feasibility of the method, which can be used in industrial production of seismic exploration.

Description

Technical Field

[0001] The present invention belongs to the field of near-surface exploration seismic data processing. Specifically, it relates to an unsupervised brittle parameter inversion method, device, and medium. Background Art

[0002] Currently, for brittle prediction of reservoirs, mainly several methods such as geophysical inversion and well logging data calculation are used. Among them, prestack AVO inversion can make full use of the geological, lithological, and fluid information contained in prestack seismic data to obtain rich elastic parameters of underground media and is widely used in brittle prediction. The theoretical basis of prestack AVO inversion is the Zoeppritz equation. However, due to the relatively complex form of this equation and inconvenience in utilization, different scholars have made approximate treatments on it. Aki and Richards expressed the reflection coefficient as a relational expression of longitudinal wave velocity, transverse wave velocity, and density, providing a theoretical basis for subsequent AVO inversion. Goodway proposed a new formula for expressing brittleness using Lame parameters and shear modulus. Zong Zhaoyun deduced the YPD approximate equation based on the Aki&Richard approximation formula, established a linear relationship between the longitudinal wave reflection coefficient and the Young's modulus reflection coefficient, Poisson's ratio reflection coefficient, and density reflection coefficient, and provided a method for directly inverting Young's modulus and Poisson's ratio. Guo et al. proposed a new formula for expressing brittleness using Young's modulus and Poisson's ratio. Zhang Guangzhi deduced approximate equations for longitudinal wave and converted wave reflection coefficients based on Eρ, Poisson's ratio, and density on the basis of the YPD equation, providing a new method for direct inversion. Since seismic inversion is a typical ill-posed problem with ill-conditioning, many scholars have carried out research on inversion algorithms, such as L1 regularization, L2 regularization methods, etc. These traditional algorithms can improve the inversion accuracy to a certain extent, but they are time-consuming in practical applications. Seismic prestack inversion needs to continue to explore in new directions. Summary of the Invention

[0003] The object of the present invention is to overcome the technical problems proposed in the background art and provide an unsupervised brittle parameter inversion method, device, and medium to improve the accuracy and efficiency of traditional prestack brittle inversion.

[0004] The specific technical solution of the present invention is as follows:

[0005] According to the first aspect of the present invention, an unsupervised brittle parameter inversion method is provided. The method includes:

[0006] Obtain prestack seismic data, and perform partial stacking to obtain multi-angle seismic data S, and send the multi-angle seismic data S into a deep learning network for training to obtain predicted brittle parameter P L ;

[0007] Construct a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; use the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D to establish a linear forward modeling equation operator A = WBD;

[0008] Multiply the predicted brittleness parameter P L by the forward modeling equation operator A to obtain the forward modeled seismic data S`, calculate the error between the forward modeled seismic data S` and the original seismic data S and backpropagate it to the deep learning network. At the same time, backpropagate a part of the error between the predicted brittleness parameter P L and the low-frequency initial model L0 to the deep learning network for constraint;

[0009] Train the network according to the backpropagated error between the forward modeled seismic data S` and the original seismic data S, iteratively update the network parameters w, and make the error between S and S` reach a preset threshold to obtain the optimal network parameters; use the optimal network parameters to optimize the deep learning network and directly predict the brittleness parameter based on the optimized deep learning network.

[0010] Furthermore, the deep learning network is a Fastformer network.

[0011] Furthermore, the step of constructing a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; using the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D to establish a linear forward modeling equation operator A = WBD specifically includes:

[0012] Based on the reflection coefficient approximation equation of E ρ , Poisson's ratio, and density:

[0013] R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1)

[0014] where θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the transverse wave velocity, R(θ) is the seismic reflection coefficient, R_Eρ is the Eρ reflection coefficient, R_δ is the Poisson's ratio reflection coefficient, R_ρ is the density reflection coefficient;

[0015] According to Equation (1), the reflection coefficient of Poisson's ratio δ is approximately expressed as:

[0016]

[0017] where i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative;

[0018] Denote L i as δi The natural logarithm of is represented by the matrix as \(R_{\delta}=DL\) δ , where: \(D\) - difference matrix; \(L\) δ — the natural logarithm of the Poisson ratio; using the angle-dependent wavelet matrix \(W(\theta)\), then the angle-dependent forward seismic record \(S(\theta)\) is expressed as:

[0019]

[0020] In the formula, \(L\) Eρ , \(L\) δ and \(L\) ρ respectively represent the natural logarithm of density volume, Poisson ratio, and density; the matrices \(a(\theta\) i ), \(b(\theta\) i ), and \(c(\theta\) i ) are all diagonal matrices. Defining the diagonal matrix as \(B\), the forward equation based on the linear approximation formula is abbreviated as:

[0021] \(S = WBDL = AL, A = WBD\ (4)\)

[0022] Furthermore, the loss function of the deep learning network is:

[0023] \(Loss = ||AP\) L - S||_2 + \mu * ||P L - L_0||_2\ (5)\)

[0024] In the formula, \(|| ||_2\) is the \(L_2\) norm, and \(\mu\) is the low-frequency constraint term parameter, taking values from 0.001 to 0.0001.

[0025] Furthermore, based on the loss function shown in formula (5), adding \(L_1\) regularization, the final loss function of the deep learning network is expressed as:

[0026] \(Loss = ||AP\) L - S||_2 + \mu * ||P L - L_0||_2 + \lambda * ||w||_1\ (6)\)

[0027] In the formula, \(\lambda\) is the regularization parameter, which decreases as the number of training rounds increases, and \(w\) is the network parameter.

[0028] According to the second aspect of the present invention, an unsupervised brittle parameter inversion device is provided, and the device includes:

[0029] A training unit, configured to obtain pre-stack seismic data, perform partial stacking to obtain multi-angle seismic data \(S\), and send the multi-angle seismic data \(S\) into a deep learning network for training to obtain predicted brittle parameters \(P\)L ;

[0030] An operator building unit, configured to construct a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; use the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D to establish a linear forward modeling equation operator A = WBD;

[0031] An error determination unit, configured to multiply the predicted brittleness parameter P L by the forward modeling equation operator A to obtain a forward modeled seismic data S`, calculate the error between the forward modeled seismic data S` and the original seismic data S and feed it back to the deep learning network, and at the same time feed back part of the error between the predicted brittleness parameter P L and the low-frequency initial model L0 to the deep learning network for constraint;

[0032] A prediction unit, configured to train the network according to the error between the fed-back forward modeled seismic data S` and the original seismic data S, iteratively update the network parameters w, make the error between S and S` reach a preset threshold, and obtain the optimal network parameters; use the optimal network parameters to optimize the deep learning network and directly predict the brittleness parameter based on the optimized deep learning network.

[0033] Further, the deep learning network is a Fastformer network.

[0034] Further, the operator building unit is further configured to:

[0035] Based on the reflection coefficient approximation equation of Eρ, Poisson's ratio, and density:

[0036] R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1)

[0037] where θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the transverse wave velocity, R(θ) is the seismic reflection coefficient, R_Eρ is the Eρ reflection coefficient, R_δ is the Poisson's ratio reflection coefficient, R_ρ is the density reflection coefficient;

[0038] According to Equation (1), the reflection coefficient approximation of Poisson's ratio δ is expressed as:

[0039]

[0040] where i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative;

[0041] Denote L i as the natural logarithm of δ i and there is It is represented by a matrix as \(R_{\delta}=DL\) δ , where: \(D\) - difference matrix; \(L\) δ — natural logarithm of Poisson's ratio; Using the angle - related wavelet matrix \(W(\theta)\), the angle - related forward seismic record \(S(\theta)\) is expressed as:

[0042]

[0043] In the formula, \(L\) Eρ , \(L\) δ and \(L\) ρ represent the density volume, Poisson's ratio, and natural logarithm of density respectively; The matrices \(a(\theta\) i ), \(b(\theta\) i ) and \(c(\theta\) i ) are all diagonal matrices. Defining the diagonal matrix as \(B\), the forward equation based on the linear approximation formula is abbreviated as:

[0044] \(S = WBDL=AL\), \(A = WBD\ (4)\)

[0045] Furthermore, the prediction unit is further configured to update the network parameter \(w\) through the following loss function:

[0046] \(Loss=\|AP\) L - S\|^2+\mu * \|P L - L_0\|^2+\lambda*\|w\|_1\ (6)

[0047] In the formula, \(\lambda\) is the regularization parameter, which decreases as the number of training rounds increases, and \(w\) is the network parameter.

[0048] According to the third aspect of the present invention, there is provided a computer - readable storage medium, on which computer - readable instructions are stored. When the computer - readable instructions are executed by a processor of a computer, the computer is caused to execute the unsupervised brittle parameter inversion method as described in each embodiment of the present invention.

[0049] According to the unsupervised brittle parameter inversion method, device, and medium provided by each embodiment of the present invention, it has at least the following beneficial effects:

[0050] First, the reliability and advancement of the method principle. It mainly uses the approximate brittle parameter inversion equation with a relatively mature theoretical basis and the Fastformer network with good application effects. It fully utilizes the ability of the deep - learning network to mine data features, drives the network prediction results with the forward equation and constrains them with the low - frequency model, combines the deep - learning framework with the physical equation to improve the prediction accuracy and efficiency of the brittle parameter inversion task, and has high reliability and advancement in principle. Compared with the conventional inversion algorithm, the prediction accuracy is higher.

[0051] Second, it does not require training labels and has strong applicability. During the entire training process, instead of using real brittleness parameters as the training set to train the network, the original seismic data is used as the input and actual labels of the network, achieving an approximate unsupervised learning effect. It is applicable to approximate pre-stack parameter inversion and actual data inversion and has great application potential.

[0052] Third, it has high computational efficiency. The inversion network selects the Fastformer network with a lightweight structure, which has low time and space complexity, fast training and convergence speed. At the same time, by making full use of GPU parallelism and CUDA acceleration, large-scale data can be processed, and it has high computing efficiency. Description of the Drawings

[0053] In order to more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the following will briefly introduce the drawings required for the description of the specific embodiments or the prior art. In all the drawings, similar elements or parts are generally identified by similar reference numerals. In the drawings, the elements or parts are not necessarily drawn to scale.

[0054] Figure 1 It is a processing flow chart of an unsupervised seismic brittleness parameter prediction method according to an embodiment of the present invention.

[0055] Figure 2 It is a multi-angle theoretical seismic record diagram, where (a)-(c) respectively represent the 10°, 20°, and 30° multi-angle theoretical seismic records.

[0056] Figure 3 It is a schematic diagram of the theoretical parameters to be inverted after taking the logarithm, where (a)-(c) respectively represent the theoretical Eρ, δ, and ρ to be inverted after taking the logarithm.

[0057] Figure 4 It is a schematic diagram of the low-frequency initial parameters, where (a)-(c) respectively represent the low-frequency initial Eρ, δ, and ρ.

[0058] Figure 5 It is the result of the inversion of the theoretical data by the present invention, where (a)-(c) respectively represent the Eρ, δ, and ρ inverted from the theoretical data by the present invention.

[0059] Figure 6 It is the result of the inversion of the theoretical data by the traditional L1 method, where (a)-(c) respectively represent the Eρ, δ, and ρ inverted from the theoretical data by the traditional L1 method.

[0060] Figure 7is the absolute error between the inversion result and the theoretical value to be inverted, where (a)-(c) respectively represent the absolute errors of Eρ, δ, ρ inverted by the present invention from the theoretical data and the theoretical value to be inverted, and (d)-(f) respectively represent the absolute errors of Eρ, δ, ρ inverted by the traditional L1 method from the theoretical data and the theoretical value to be inverted.

[0061] Figure 8 is the seismic record diagram of the actual work area from multiple angles, where (a)-(c) respectively represent the seismic records of the actual work area from multiple angles of 12°, 24°, and 36°.

[0062] Figure 9 is the schematic diagram of the low-frequency initial parameters interpolated by wells in the actual work area, where (a)-(c) respectively represent the low-frequency initial Eρ, δ, ρ interpolated by wells in the actual work area.

[0063] Figure 10 is the parameter effect diagram directly inverted by the present invention in the actual work area, where (a)-(c) respectively represent Eρ, δ, ρ.

[0064] Figure 11 is the parameter effect diagram directly inverted by the traditional L1 method in the actual work area, where (a)-(c) respectively represent Eρ, δ, ρ.

[0065] Figure 12 is the parameter effect diagram indirectly inverted by the software in the actual work area, where (a)-(c) respectively represent Eρ, δ, ρ.

[0066] Figure 13 is the training flow chart of the unsupervised method driven by physical equations of the present invention.

[0067] Figure 14 is the structure diagram of an unsupervised seismic brittleness parameter prediction device according to an embodiment of the present invention. Detailed implementation manners

[0068] Next, the technical solutions in the embodiments of the present invention will be described clearly and completely. Apparently, the described embodiments are only a part of the embodiments of the present invention, rather than all of them. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.

[0069] In order to make the purpose, technical solutions and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention, rather than to limit the invention.

[0070] Now, the present invention will be further described with reference to the accompanying drawings of the specification.

[0071] At present, various traditional pre-stack seismic inversion methods have certain limitations in use. In recent years, machine learning and deep learning have shown great potential when applied to various geophysical research problems, providing automated performance in some tasks. To solve the problems of low efficiency and ill-posedness in traditional pre-stack brittle inversion, deep learning (DL) is introduced into pre-stack seismic brittle inversion. Since there is very little labeled data for network training in actual data, the dependence on a large number of training samples affects the application of the (DL) method in pre-stack seismic inversion. Therefore, the embodiment of the present invention proposes an unsupervised brittle parameter inversion method, specifically an unsupervised DL pre-stack seismic brittle parameter inversion method driven by physical equations. This method integrates the DL framework and an inversion method based on the brittle parameter inversion equation, providing a solution strategy that does not require training labels for DL pre-stack seismic inversion work, while improving the inversion accuracy and efficiency, and is applicable to large-scale data processing.

[0072] The development of shale reservoirs mainly depends on hydraulic fracturing, where brittleness is an important parameter characterizing shale quality and fracturability. Pre-stack seismic inversion can accurately obtain the brittle parameters of underground reservoir media, and the method proposed in the embodiment of the present invention can improve the accuracy and efficiency of traditional pre-stack brittle inversion.

[0073] As Figure 13 shown, the basic technical principle of this method is as follows:

[0074] Conventional supervised deep learning takes seismic data as the input of the neural network and the corresponding seismic inversion parameters as labels. A large number of established training data sets are used to train the neural network, and then the trained network model is used to predict the seismic data to achieve the purpose of inversion. The present invention integrates the DL framework and the brittle parameter inversion equation, providing a solution strategy that does not require training labels for DL pre-stack seismic inversion work. First, the original seismic data is input into the network for brittle parameter prediction, and a low-frequency model is introduced as a physical constraint. Then, an approximate unsupervised inversion method is proposed. This method inserts a forward modeling equation of linear approximation as a module into the network framework, and uses the predicted data of the network as the input of the forward modeling module to obtain the corresponding forward modeled seismic data. Through the backpropagation of the error between the forward modeled seismic data and the original seismic data and the error between the network predicted data and the low-frequency model, the network parameters are iteratively optimized to minimize the error, and the network is used to directly predict the brittle parameters. During the entire training process, instead of using real brittle parameters as the training set to train the network, the original seismic data is used as the input and actual label of the network, achieving the effect of approximate unsupervised learning.

[0075] Based on the above technical principle, please refer to Figure 1, an unsupervised brittle parameter inversion method proposed in this embodiment. This method starts from step S100, obtains pre-stack seismic data, and performs partial stacking to obtain multi-angle seismic data S. The multi-angle seismic data S is sent into a deep learning network for training to obtain a predicted brittle parameter P L ;

[0076] In some embodiments, the deep learning network is a Fastformer network. Compared with recurrent neural network (RNN, LSTM, GRU, etc.) networks, FastFormer is a variant of Transformer that can achieve context modeling with linear complexity. Its time and space complexity are more efficient than that of the standard Transformer, and the total number of training parameters is also greatly reduced. The efficiency of pre-stack seismic inversion for long-time series seismic data is relatively low. Therefore, the Fastformer network is more suitable for our task than conventional recurrent neural networks and other Transformer variant networks

[0077] In step S200, according to the multi-angle seismic data S, a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D are constructed; the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D are used to establish a linear forward equation operator A = WBD

[0078] In some embodiments, the constructing a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; using the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D to establish a linear forward equation operator A = WBD specifically includes:

[0079] Based on the reflection coefficient approximation equation of Eρ, Poisson's ratio, and density:

[0080] R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1)

[0081] In the formula, θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the transverse wave velocity, R(θ) is the seismic reflection coefficient R_Eρ is the Eρ reflection coefficient R_δ is the Poisson's ratio reflection coefficient R_ρ is the density reflection coefficient

[0082] According to formula (1), the reflection coefficient of Poisson's ratio δ is approximately expressed as:

[0083]

[0084] In the formula, i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative

[0085] Let L i be the natural logarithm of δ i , then we have It is represented by a matrix as R_δ = DL δ , where: D - difference matrix; L δ —the natural logarithm of the Poisson ratio; using the angle - related wavelet matrix W(θ), then the angle - related forward seismic record S(θ) is expressed as:

[0086]

[0087] In the formula, L Eρ , L δ and L ρ represent the natural logarithm of density volume, Poisson ratio and density respectively; the matrices a(θ i ), b(θ i ) and c(θ i ) are all diagonal matrices. Defining the diagonal matrix as B, then the forward equation based on the linear approximation formula is abbreviated as:

[0088] S = WBDL = AL, A = WBD (4)

[0089] In step S300, multiply the predicted brittle parameter P L by the forward equation operator A to obtain the forward seismic data S`, calculate the error between the forward seismic data S` and the original seismic data S and back - propagate it to the deep learning network. At the same time, part of the error between the predicted brittle parameter P L and the low - frequency initial model L0 is back - propagated to the deep learning network for constraint.

[0090] Finally, in step S400, train the network according to the error between the back - propagated forward seismic data S` and the original seismic data S, iteratively update the network parameter w, make the error between S and S` reach the preset threshold to obtain the optimal network parameter; optimize the deep learning network using the optimal network parameter, and directly predict the brittle parameter based on the optimized deep learning network.

[0091] In some embodiments, the loss function of the deep learning network is:

[0092] Loss = ||AP L - S||2 + u * ||P L - L0||2 (5)

[0093] In the formula, || ||2 is the L2 norm, μ is the low - frequency constraint term parameter, taking values from 0.001 to 0.0001.

[0094] In some embodiments, considering that the seismic brittleness parameters are sparse, in order to improve the overall inversion effect, the present invention adds L1 regularization to the loss function to make the inversion result closer to the true value. The final loss function is:

[0095] Loss=||AP L -S||2+u * ||P L -L0||2+λ * ||w||1 (6)

[0096] In the formula, λ is the regularization parameter, which decreases as the number of training rounds increases, and w is the network parameter.

[0097] Exemplarily, the architecture of the deep learning network (Fastformer network) is built on the PyTorch platform and uses GPU parallelism to accelerate training. The network is trained using the ADAM optimizer with its hyperparameter settings (initial learning rate: LR = 1e-3, batch size: batch_size = 32, decay rate: decay_weight = 1e-5).

[0098] To verify the effect of the unsupervised seismic brittleness parameter prediction method driven by physical equations, the following analyzes are carried out by taking the inversion results of the theoretical two-dimensional model and the actual data of a certain shale work area in Sichuan as examples respectively.

[0099] To verify the inversion effect of the method, a part of the Marmousi2 theoretical model is selected to generate reflection coefficients at 3 incident angles, and the incident angles are 10°, 20°, and 30° respectively. Then, a multi-angle seismic data volume is synthesized by convolving a zero-phase wavelet with a main frequency of 30 Hz and the reflection coefficients. There are a total of 750 traces, with a sampling interval of 1 ms, and each trace includes 500 sampling points. Figure 2 (a)-(c) are the angle seismic data volumes at incident angles of 10°, 20°, and 30° respectively. Among them, in order to meet the inversion theory, the true model is logarithmically taken to obtain the theoretical model to be inverted as shown in Figure 3 and a low-frequency initial model as shown in Figure 4 is obtained on the basis of the model to be inverted. The multi-angle seismic data and the low-frequency model are fed into the network and trained for 50 Epoch, with a total of 1200 iterations. After the training is completed, the saved network model is directly used to predict the brittleness parameters.

[0100] As shown in Figure 5 , the inversion result of the present invention basically reflects the change trend of the brittleness parameters. Compared with Figure 3Comparing (a)-(b), the Eρ and Poisson's ratio directly inverted through the network have higher accuracy. The inversion results in the fault area of the model are also relatively stable, with a high degree of agreement with the true model, good continuity in the horizontal direction, and clear formation boundaries in the vertical direction. Figure 5 (c) shows the density inversion result. As can be seen from the comparison with Figure 3 (c), the density inversion result can better reflect the geological characteristics of the model. To sum up, the inversion result of the method proposed in the present invention has a smaller error and higher accuracy. To further compare the inversion effects, the traditional L1 method was used for inversion and a comparative analysis was carried out with the network inversion effect. Figure 6 (a)-(c) are the inversion results of the traditional L1 method. To more clearly show the differences between the two methods, the absolute errors between the two different inversion results and the true model to be inverted were calculated respectively. From Figure 7 the inversion error diagrams of (a)-(f), it can be observed that the inversion results of the unsupervised deep learning method have a phenomenon of reduced error in most areas. Especially in the part of the red box, it can be seen that the inversion result error of the method in this paper is significantly reduced, the result is more stable, and it is more consistent with the original model. At the same time, under the condition of the same amount of inversion data, as shown in Table 1, the present invention has higher computing efficiency than the traditional algorithm. If there is a GPU graphics card with stronger computing power, the acceleration efficiency can be further improved by using CUDA.

[0101] Table 1 Comparison of theoretical data calculation efficiency

[0102]

[0103] To further verify the effectiveness of the present invention, it was applied to the inversion of brittle parameters of actual data. The data is from a shale work area in the Sichuan Basin, China, and part of the angle stack profile is shown. Figure 8 (a)-(c) are seismic data volumes with incident angles of 12°, 24°, and 36° respectively. This profile has 800 traces of data, the sampling time period is 1200 - 1920 ms, the longitudinal time length is 720 ms, and the sampling interval is 4 ms. The upper and lower boundary lines of the profile represent two geological interpretation horizons, and the black vertical lines are the well positions.

[0104] Applying the inversion method to the above actual data, first extract the statistical wavelet from the original seismic data to construct the forward operator A, and then interpolate the multi-well data in the work area to obtain the low-frequency model such as Figure 9 (a)-(c). Send the multi-angle seismic data and the original low-frequency model after taking the logarithm into the network, train for 50 Epochs, and iterate 1200 times in total. After the training is completed, directly predict the brittle parameters using the saved network model. Figure 10 、 Figure 11 、 Figure 12The final inversion results of the present invention, the traditional L1 method, and commercial software are respectively shown. Among them, the inversion resolution of the commercial software is slightly lower, and the overall inversion result is relatively consistent with the other two methods, while the inversion results in some areas are slightly different from the other two methods. Compared with the traditional L1 method, the overall resolution of the inversion result of the unsupervised deep learning method has been well improved. For the convenience of observation, some data within the black square are magnified, and it can be clearly seen that the lateral continuity of the inversion result of the proposed method is better, and at the same time, the vertical resolution has been effectively improved, which can effectively depict the formation characteristics. Similarly, under the condition of the same amount of inversion data, as shown in Table 2, the present invention has higher operation efficiency than the traditional algorithm, which shows the advantages of the method.

[0105] Table 2 Comparison of calculation efficiency of theoretical data

[0106]

[0107] The embodiment of the present invention also provides an unsupervised brittle parameter inversion device, as Figure 14 shown, the device 1400 includes:

[0108] A training unit 1401, configured to obtain pre-stack seismic data, perform partial stacking to obtain multi-angle seismic data S, and send the multi-angle seismic data S into a deep learning network for training to obtain a predicted brittle parameter P L ;

[0109] An operator establishment unit 1402, configured to construct a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; use the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D to establish a linear forward modeling equation operator A = WBD;

[0110] An error determination unit 1403, configured to multiply the predicted brittle parameter P L by the forward modeling equation operator A to obtain a forward modeled seismic data S`, calculate the error between the forward modeled seismic data S` and the original seismic data S and send it back to the deep learning network, and at the same time send back part of the error between the predicted brittle parameter P L and the low-frequency initial model L0 to the deep learning network for constraint;

[0111] A prediction unit 1404, configured to train the network according to the error between the back-propagated forward modeled seismic data S` and the original seismic data S to iteratively update the network parameter w, make the error between S and S` reach a preset threshold to obtain the optimal network parameter; optimize the deep learning network with the optimal network parameter, and directly predict the brittle parameter based on the optimized deep learning network.

[0112] In some embodiments, the deep learning network is a Fastformer network.

[0113] In some embodiments, the operator building unit is further configured to:

[0114] Based on the reflection coefficient approximation equation of Eρ, Poisson's ratio, and density:

[0115] R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1)

[0116] Where θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the transverse wave velocity, and R(θ) is the seismic reflection coefficient. R_Eρ is the Eρ reflection coefficient. R_δ is the Poisson's ratio reflection coefficient. R_ρ is the density reflection coefficient.

[0117] According to Equation (1), the reflection coefficient of Poisson's ratio δ is approximately expressed as:

[0118]

[0119] Where i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative.

[0120] Denote L i as the natural logarithm of δ i , then there is Represented by a matrix as R_δ = DL δ , where: D - difference matrix; L δ - the natural logarithm of Poisson's ratio; Using the angle-dependent wavelet matrix W(θ), the angle-dependent forward seismic record S(θ) is expressed as:

[0121]

[0122] Where L Eρ , L δ and L ρ represent the natural logarithms of density volume, Poisson's ratio, and density respectively; The matrices a(θ i ), b(θ i ), and c(θ i ) are all diagonal matrices. Define the diagonal matrix as B, then the forward equation based on the linear approximation formula is abbreviated as:

[0123] S = WBDL = AL, A = WBD (4)

[0124] In some embodiments, the prediction unit is further configured to update the network parameter w through the following loss function:

[0125] Loss = ||APL -S||2+u*||P L -L0||2+λ*||w||1 (6)

[0126] Wherein, λ is the regularization parameter, which decreases as the number of training rounds increases during setting, and w is the network parameter.

[0127] It should be noted that the unsupervised brittle parameter inversion device provided in the embodiments of the present invention and the previously described unsupervised brittle parameter inversion method belong to the same technical concept and have the same beneficial effects, which will not be elaborated here.

[0128] The embodiments of the present invention provide a computer-readable storage medium, on which computer-readable instructions are stored. When the computer-readable instructions are executed by a processor of the computer, the computer is enabled to execute the unsupervised brittle parameter inversion method as described in various embodiments of the present invention.

[0129] The above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it; although the present invention has been described in detail with reference to the foregoing embodiments, those of ordinary skill in the art should understand that: they can still modify the technical solutions recorded in the foregoing embodiments, or perform equivalent replacements on some or all of the technical features; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the scope of the technical solutions of the embodiments of the present invention, and they should all be covered by the scope of the claims and the description of the present invention.

Claims

1. An unsupervised brittle parameter inversion method, characterized in that, The method includes: Obtain pre-stack seismic data and perform partial stacking to obtain multi-angle seismic data S, and send the multi-angle seismic data S into a deep learning network for training to obtain a predicted brittleness parameter P L ; Constructing a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; establishing a linear forward modeling equation operator A = WBD using the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D; Multiply the predicted brittleness parameter P L by the forward modeling equation operator A to obtain the forward modeled seismic data S`, calculate the error between the forward modeled seismic data S` and the original seismic data S, and backpropagate the error to the deep learning network. At the same time, part of the error between the predicted brittleness parameter P L and the low-frequency initial model L0 is backpropagated to the deep learning network for constraint; Training the network based on the error between the back-propagated forward modeling seismic data S' and the original seismic data S to iteratively update the network parameter w, making the error between S and S' reach a preset threshold to obtain the optimal network parameter; optimizing the deep learning network using the optimal network parameter, and directly predicting the brittleness parameter based on the optimized deep learning network.

2. The method according to claim 1, characterized in that, The deep learning network is a Fastformer network.

3. The method according to claim 1, characterized in that, The constructing a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; establishing a linear forward modeling equation operator A = WBD using the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D specifically includes: Based on the reflection coefficient approximation equation of density volume Eρ, Poisson's ratio, and density: R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1) where θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the transverse wave velocity, and R(θ) is the seismic reflection coefficient, R_Eρ is the Eρ reflection coefficient, R_δ is the Poisson's ratio reflection coefficient, R_ρ is the density reflection coefficient; According to Equation (1), the reflection coefficient of Poisson's ratio δ is approximately expressed as: where i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative; Denote L i as δ i the natural logarithm of, there is represented by a matrix as R_δ = DL δ , where: D - difference matrix; L δ — the natural logarithm of the Poisson ratio; using the angle-related wavelet matrix W(θ), then the angle-related forward seismic record S(θ) is expressed as: Wherein, L Eρ , L δ and L ρ represent the density volume, Poisson's ratio, and the natural logarithm of density, respectively; the matrices a(θ i ), b(θ i ), and c(θ i ) are all diagonal matrices. Defining the diagonal matrix as B, the forward equation based on the linear approximation formula is abbreviated as: S = WBDL = AL, A = WBD (4).

4. The method according to claim 1, characterized in that, The loss function of the deep learning network is: Loss=||AP L -S||2+u*||P L -L0||2 (5) where ||||2 is the L2 norm, μ is the parameter of the low-frequency constraint term, taking 0.001 - 0.0001.

5. The method according to claim 4, characterized in that, Based on the loss function shown in Equation (5), adding L1 regularization, the final loss function of the deep learning network is expressed as: Loss=||AP L -S||2+u*||P L -L0||2+λ*||w||1 (6) where λ is the regularization parameter, which decreases as the number of training rounds increases, and w is the network parameter.

6. An unsupervised brittle parameter inversion device, characterized in that, The device includes: A training unit, configured to obtain pre-stack seismic data, perform partial stacking to obtain multi-angle seismic data S, and send the multi-angle seismic data S into a deep learning network for training to obtain a predicted brittleness parameter P L ; An operator establishment unit configured to construct a low-frequency initial model L0, a wavelet kernel matrix W, an angle coefficient matrix B, and a difference matrix D according to the multi-angle seismic data S; establish a linear forward modeling equation operator A = WBD using the wavelet kernel matrix W, the angle coefficient matrix B, and the difference matrix D; An error determination unit configured to perform matrix multiplication on the predicted brittleness parameter P L and the forward equation operator A to obtain forward seismic data S`, calculate the error between the forward seismic data S` and the original seismic data S, and feed the error back to the deep learning network. At the same time, the predicted brittleness parameter P L and a partial error of the low-frequency initial model L0 are fed back to the deep learning network for constraint; A prediction unit configured to train the network based on the error between the back-propagated forward modeling seismic data S' and the original seismic data S to iteratively update the network parameter w, making the error between S and S' reach a preset threshold to obtain the optimal network parameter; optimize the deep learning network using the optimal network parameter, and directly predict the brittleness parameter based on the optimized deep learning network.

7. The device according to claim 6, characterized in that, The deep learning network is a Fastformer network.

8. The device according to claim 6, characterized in that, The operator establishment unit is further configured to: Based on the reflection coefficient approximation equation of Eρ, Poisson's ratio, and density: R(θ) = a(θ)R_Eρ + b(θ)R_σ + c(θ)R_ρ (1) where θ is the azimuth angle, k is the ratio of the longitudinal wave velocity to the shear wave velocity, and R(θ) is the seismic reflection coefficient. R_Eρ is the Eρ reflection coefficient. R_δ is the Poisson's ratio reflection coefficient. R_ρ is the density reflection coefficient. According to Equation (1), the reflection coefficient of Poisson's ratio δ is approximately expressed as: where i is 1, 2, 3..., δ is Poisson's ratio, and Δ is the derivative; Denote L i as the natural logarithm of δ i , then we have It is represented by the matrix as \(R_{\delta}=DL\) δ , where: D - difference matrix; L δ — the natural logarithm of the Poisson ratio; Using the angle-dependent wavelet matrix \(W( heta)\), then the angle-dependent forward seismic record \(S( heta)\) is expressed as: where L Eρ , L δ and L ρ represent the density volume, Poisson's ratio, and natural logarithm of density, respectively; the matrices a(θ i ), b(θ i ), and c(θ i ) are all diagonal matrices. Defining the diagonal matrix as B, the forward equation based on the linear approximation formula is abbreviated as: S = WBDL = AL, A = WBD (4).

9. The device according to claim 6, wherein, The prediction unit is further configured to update the network parameter w through the following loss function: Loss=||AP L -S||2+u*||P L -L0||2+λ*||w||1 (6) Wherein, λ is a regularization parameter, which decreases as the number of training rounds increases during setting, and w is a network parameter.

10. A computer-readable storage medium, wherein, It stores computer-readable instructions, which, when executed by a processor of a computer, cause the computer to execute the method according to any one of claims 1-5.

Images