A VSP acoustic full waveform inversion method of a double-branch physically driven recurrent neural network

Abstract

The application discloses a VSP acoustic wave full waveform inversion method of a double-branch physical driving cyclic neural network, applied to the field of seismic data processing, and aims at the problems of weak explanation and high calculation cost of the existing seismic velocity inversion technology in a full data driving mode. Since a wave equation describes the change of a wave field in time and space, a time stepping method is adopted in numerical simulation, that is, the wave field state of the next moment is calculated according to the state of the previous moment. Thus, the wave field state of the previous moment is taken as a hidden layer of the RNN, the acoustic wave equation is solved by finite difference, and the velocity parameter in the equation is taken as a trainable parameter in the physical driving RNN forward network. The shot record output by the physical driving RNN forward network each time is taken as a predicted value. The shot record output by the physical driving RNN forward network by setting a real velocity parameter is taken as an observed value. The process of correcting the velocity parameter by reverse propagation of the loss between the predicted value and the observed value is called a velocity inversion process.

Description

Technical Field

[0001] The present invention belongs to the field of seismic data processing, and in particular relates to a full waveform inversion technology. Background Art

[0002] Full Waveform Inversion (FWI) uses the maximum wavelength information contained in seismic data to reconstruct underground parameters by minimizing the error between observed and predicted data. Full waveform inversion was first proposed by Lailly and Tarantola. The full waveform inversion problem is a strongly nonlinear problem, so FWI is highly dependent on the initial model. In addition, the data acquisition during seismic exploration lacks low-frequency components, which makes the inversion prone to falling into local minima and cycle jump problems. Warner and Guasch proposed adaptive waveform inversion (AWI) to reduce dependence on data and cycle jump problems, but it may reduce the resolution of the inversion results and increase computational complexity. In general, the FWI method faces challenges such as how to balance computational efficiency, result accuracy and dependence on the initial model.

[0003] Full waveform inversion is a strongly nonlinear method, and its results are multi-knotted. Therefore, it is highly dependent on the accuracy of the initial data. The key to traditional full waveform inversion is how to optimize the objective function, that is, the error between the observed data and the predicted data. Generally, the adjoint state method is used to calculate the gradient of the objective function with respect to the velocity parameter, and then a suitable optimization algorithm is selected to update the velocity parameter. Common optimization algorithms include the steepest descent method, conjugate gradient method, Newton method, trust region method, etc.

[0004] With the advancement of large-scale computing power and the rapid development of deep learning technology, deep learning has become increasingly widely used in geophysical inversion. Deep learning has demonstrated strong potential in addressing nonlinear problems, providing a novel research direction for seismic velocity inversion. Traditional deep learning inversion methods generally adopt a fully data-driven approach, using a pre-constructed training dataset to allow a neural network to learn the mapping between seismic observations and geophysical models. This approach makes deep learning easy to implement in seismic exploration, but it also faces challenges such as limited interpretability and computational cost. In recent years, physics-inspired neural networks (PINNs) have begun to gain widespread application. These methods integrate the forward module of wave equation physics modeling into the network loop and use the physics-driven data residual as the loss function for full waveform inversion. Compared to fully data-driven deep learning methods, this approach is more flexible, can incorporate structural information from multiple data types and any structure, and has stronger noise immunity, thereby improving the accuracy and reliability of inversion results. However, PIRNN full waveform inversion often consumes significant GPU computing resources, especially when dealing with large-scale three-dimensional problems. SUMMARY

[0005] To solve the above technical problems, the application provides a VSP acoustic full waveform inversion method of a dual-branch physical driving recurrent neural network, which combines a physical driving RNN full waveform inversion network with scale decomposition, establishes a Dual-PIRNN (Dual Physical-Informed Recurrent Neural Network based on a dual-branch physical driving recurrent neural network), and realizes a low-GPU resource consumption VSP acoustic full waveform inversion method based on machine learning and physical guidance.

[0006] The technical scheme adopted by the application is as follows: a VSP acoustic full waveform inversion method of a dual-branch physical driving recurrent neural network, which constructs a dual-branch physical driving recurrent neural network, and the dual-branch physical driving recurrent neural network specifically comprises: an original scale real velocity model forward branch, an original scale initial velocity model forward branch, and a reduced scale velocity model inversion branch; the original scale real velocity model forward branch, the original scale initial velocity model forward branch, and the reduced scale velocity model inversion branch are realized by using a physical driving RNN forward network.

[0007] The acoustic wave equation is solved by finite difference, and the velocity parameter in the equation is used as a trainable parameter in the physical driving RNN forward network; the output of the original scale real velocity model forward branch is used as an observation value, the output of the original scale initial velocity model forward branch is used as a prediction value, the loss between the prediction value and the observation value is input into the reduced scale velocity model inversion branch, the velocity parameter is updated, and the velocity parameter on the reduced scale velocity model inversion branch is updated to the original scale initial velocity model forward branch by upsampling at the end of each iteration inversion to perform forward again and update cyclically.

[0008] The application has the following beneficial effects: the Dual-PIRNN full waveform inversion network structure proposed in the application combines physical wave equation constraints, ensures the consistency of the synthesized seismic data and the observed seismic data, and thus improves the accuracy of the inversion result. The application introduces scale decomposition and constructs a dual-branch physical driving recurrent neural network full waveform inversion architecture, effectively reduces the GPU resource consumption in inversion, and at the same time ensures the accuracy of the inversion result. The scale decomposition in the application can reduce local minima in inversion and reduce the situation of falling into periodic jumps, so that the velocity model can be obtained faster. BRIEF DESCRIPTION OF DRAWINGS

[0009] Figure 1 It is a spread form of the PIRNN architecture.

[0010] Figure 2 It is an internal structure of the staggered grid finite difference RNN operator.

[0011] Figure 3 This is the Dual-PRINN full waveform inversion flow chart;

[0012] Figure 4 are the correct model and the initial model;

[0013] Among them, (a) is the true velocity model, (b) is the initial model under Gaussian smoothing with variance sigma = 5, and (c) is the initial model under Gaussian smoothing with variance sigma = 10;

[0014] Figure 5 is the inversion result;

[0015] Among them, (a) is the velocity model obtained under the initial model of sigma = 5 and the time scale inversion, (b) is the velocity model obtained under the initial model of sigma = 5 and the time and space scale inversion, (c) is the velocity model obtained under the initial model of sigma = 10 and the time scale inversion, and (d) is the velocity model obtained under the initial model of sigma = 10 and the time and space scale inversion;

[0016] Figure 6 is the velocity profile;

[0017] Among them, (a) is the velocity curve of the initial model with different smoothing degrees after time-reduced inversion at x=300m, (b) is the velocity curve of the initial model with different smoothing degrees after time-reduced inversion at x=700m, (c) is the velocity curve of the initial model with different smoothing degrees after time-reduced and space-scaled inversion at x=300m, and (d) is the velocity curve of the initial model with different smoothing degrees after time-reduced and space-scaled inversion at x=700m. DETAILED DESCRIPTION

[0018] To facilitate those skilled in the art to understand the technical content of the present invention, the present invention is further explained below with reference to the accompanying drawings.

[0019] 1. First, let’s explain the recurrent neural network:

[0020] Recurrent Neural Networks (RNN) is a type of artificial neural network. RNN can be divided into two categories: finite impulse response networks with directed acyclic graphs and infinite impulse response networks with directed cyclic graphs. Both types of networks can represent dynamic time behaviors. Unlike deep neural networks, RNN can use its internal state (memory) to process input sequences so that the current output results are affected by past results. This makes them suitable for time-dependent signal processing tasks, such as text and audio processing. Since the wave equation describes the changes in the wave field in time and space, the time-stepping method is generally used when performing numerical simulation solutions, that is, the wave field state at the next moment is calculated based on the wave field state at the previous moment. The core feature of RNN is that it can process time series data, and its hidden state can carry the information of the previous moment to the next moment, which is highly consistent with the time-stepping solution method of the wave equation. Therefore, the seismic forward modeling calculation (time-stepping solution of the wave equation) can be integrated into the RNN unit through finite differences. The RNN network forward modeling architecture is as follows: Figure 1 As shown in the figure, it takes the wave field information of the previous moment as input, generates a wave field record at the current moment, and saves the wave field record generated at this moment in the GPU and passes it to the next moment to continue generating wave field records. Under the VSP Vertical Seismic Profiling (vertical seismic profiling) observation system, the space-time field information about the velocity model is obtained by collecting the wave field records at each time point.

[0021] 2. Then, the staggered grid finite difference RNN operator network structure adopted by the present invention is described:

[0022] The solution to the first-order velocity-pressure acoustic wave equation adopts a high-order finite difference format with a staggered grid form. Existing research has shown that the process of integrating seismic forward modeling into the forward process of the RNN network to achieve full waveform inversion is equivalent to the traditional full waveform inversion method. The velocity parameters are used as trainable parameters, and seismic forward modeling is performed based on the velocity parameters. The shot collection records at each time point are collected as predicted values. The error between the observed value and the predicted value is used for back propagation and updating the velocity parameters. This process is called the velocity inversion process. The automatic differentiation mechanism of the neural network is used, and other information can be easily added under the PIRNN network to improve the inversion accuracy and robustness. The internal structure of the staggered grid finite difference RNN operator is as follows: Figure 2 shown.

[0023] Figure 2The operator demonstrated here is a high-precision forward modeling operator based on PIRNN (Physically Informed Recurrent Neural Network), with the input being velocity parameters, source location, and wavefield information from the previous moment, and the output being the VSP shot gather record at the current moment. is the input at the current moment. In this structure, As the internal parameters of each SGFD RNN Operator (Staggered Grid Finite Difference Recurrent Neural Network Operator) unit, it is used to propagate between different RNN units. It contains the chain relationship related to the speed parameter v at each moment. The relationship in the entire time field will be calculated in the backpropagation according to the automatic differentiation property with respect to v, and then the speed parameter is updated through the Adam optimizer.

[0024] This embodiment takes the first-order pressure-velocity fluctuation equation as an example to construct the RNN recursive formula with second-order time accuracy in the time domain and second-order staggered grid finite difference accuracy:

[0025]

[0026] where p x and p z are the pressure components in the x and z directions, v x ,v z is the particle velocity in the x and z directions, and v is the speed of sound, which is the desired inversion. Accurate forward modeling is the basis for full waveform inversion. This paper solves the first-order pressure-velocity acoustic wave equation using high-order finite differences in a staggered grid format, while adding SPML boundaries to ensure high-precision forward modeling.

[0027] Discretize the above equation and replace the pressure p x ,p z Translated to the integer grid points, the particle velocity v x , v z Defined on half grid points, Expanding at time and [ixΔt,izΔz], we get the following formula:

[0028]

[0029]

[0030] Where △x represents the discrete step length in the x direction of the finite difference space, △z represents the discrete step length in the z direction of the finite difference space, and △x and △z are also called grid spacing; △t represents the discrete step length in time, which represents the interval between two adjacent time points in the grid;

[0031] Further processing of the above discrete equations yields the following equations about p: x , p z ,v x , v z The recursive formula is:

[0032]

[0033] Combine Figure 2 It can be seen that in each time cycle, the present invention first uses the pressure p at time kΔt k and Particle velocity at time and Obtain Particle velocity at time and Then according to the obtained Particle velocity at time and and the pressure of kΔt and To calculate the pressure at time (k+1)Δt and In this process, each time step to solve the wave equation uses the previous moment kΔt and This process allows seismic forward modeling to be performed through RNN, where each layer of the RNN can be constructed to calculate and store the wave field information at a certain moment. This structure allows us to map the time evolution process of the wave field to the hierarchical structure of the RNN. Specifically, each layer of the RNN can represent the state of the wave field at a specific time step, including pressure and particle velocity. Under this framework, each layer of the RNN not only stores the wave field information at the current moment, but also receives the output of the previous layer as input. This structure simulates the time dependence of the wave equation well, where the solution at each moment depends on the state at the previous moment. The present invention sets the above-mentioned forward solution process as the SGFD RNNOperator operator for forward solution at each moment.

[0034] 3. Explain the dual-branch physics-driven recurrent neural network based on the physical embedding recurrent neural network:

[0035] The VSP acoustic full waveform inversion structure of the dual-branch physical driven recurrent neural network based on the physical embedded recurrent neural network is as follows: Figure 3First, there is a real original scale velocity model, and the VSP observation shot record d is obtained by forward modeling with the SGFD RNN Operator. obs In each iterative inversion process, the present invention first performs forward modeling on the original scale initial model through SGFD RNN Operator to obtain the simulated shot data d cal , through the L1 norm of d obs and d cal Find the loss and get the loss on the original scale Loss original The loss is used for inversion on the downscaled branch model through scale conversion to update the velocity parameters. After each iterative inversion, the velocity parameters on the downscaled model are updated to the original scale model through upsampling and forward modeling is performed again, and the update cycle is repeated.

[0036] Compared to completely data-driven deep learning methods, the PIRNN-based full waveform inversion method is more flexible and can incorporate structural information from multiple data types and any structure. It has stronger noise immunity, thereby improving the accuracy and reliability of the inversion results. However, the PIRNN full waveform inversion method often consumes a large amount of GPU computing resources, especially when dealing with large-scale three-dimensional problems.

[0037] The Dual-PIRNN full-waveform inversion network structure proposed in this paper incorporates physical wave equation constraints to ensure consistency between synthetic and observed seismic data, thereby improving the accuracy of the inversion results. This paper introduces scale decomposition and constructs a dual-branch, physically driven, recurrent neural network full-waveform inversion architecture, effectively reducing GPU resource consumption during inversion while ensuring the accuracy of the inversion results. The scale decomposition in this paper can reduce local minima in the inversion and reduce the risk of cycle jumps, thereby faster convergence to the velocity model.

[0038] The implementation process of the present invention includes the following contents:

[0039] 1. Staggered Grid Finite Difference Acoustic Wave Equation Embedded in a Recurrent Neural Network

[0040] In the time domain, the first-order velocity-pressure form of the acoustic wave equation is:

[0041]

[0042] Where p is the pressure, p x and p z Represents the pressure in the x and z directions, v x ,v z are the x and z components of the particle velocity, v is the speed of sound, and (r, t) represents the value of the physical quantity at position r at time t.

[0043] In order to make the network suitable for velocity model inversion, the RNN structure contains trainable parameters, such as Figure 2 As shown, that is, the velocity parameter v. Then, the seismic shot gather record (observation record) at the wave field detector of the network simulation at time t obtained according to the set correct velocity model is used as a label; the initial velocity model is set with reference to the correct velocity model, and can be smoothed or globally offset. The seismic shot gather obtained by inputting the initial velocity model into the network is used as the predicted seismic record. Information is passed through Figure 1 The trainable RNN network structure shown completes forward propagation and simulates the forward modeling of the discrete wave field of the acoustic wave model.

[0044] The standard network training problem is to minimize the error between the network output (predicted seismic records) and the labeled data (observed seismic records). The trainable network weights are discrete, representing the sampled values ​​of the velocity field. The inversion process of the Dual-PIRNN can be formulated as an optimization problem that minimizes the difference between the observed and predicted data. The L2 norm is chosen to measure the data difference:

[0045]

[0046] Where, ‖·‖ 2 represents the L2 norm, Ns is the number of earthquake sources, and T is the total time. cal and d obs Representing predicted data and observed data, respectively. In a deep learning framework, the relationships between variables are automatically recorded. During backpropagation, gradients are calculated through automatic differentiation, and parameters are updated using an optimizer, such as the Adam optimizer.

[0047] 2. Dual-PIRNN network structure

[0048] By converting the loss obtained by forward modeling under the original scale branch to the downscale branch according to the scale, the high resource consumption in the PIRNN inversion can be effectively reduced. Specifically, the difference between the original scale model and the downscale model lies in the different model parameter settings, such as the time sampling rate dt of the source function, the number of time points NT of the numerical simulation, and the multiple relationship between the number of grids. In each iterative inversion process, such as Figure 2 As shown, first, the SGFDRNN Operator is used to perform forward modeling on the original scale initial model to obtain the simulated shot gather data d cal , through the L1 norm of d obs and d cal Find the loss and get the loss on the original scale Loss original :

[0049] Loss original =dobs -d cal

[0050] Here, the present invention obtains d by the VSP observation system cal and d obs Is a time series matrix of shape NTxnz, so the calculated Loss original It is also a time series matrix with the shape of NT×nz. Then, the present invention uses the downsampling method based on average pooling to reduce the loss original Extract the loss that can match the downscaling model dowm Finally, the VSP shot gather record d on the downscaled model is obtained by forward modeling on the other downscaled model through the SGFD RNN Operator. Dcal . down and d Dcal Add together to get the observation record d on the downscaled size model Dobs :

[0051] Loss down =ExtractionOperator(Loss original )

[0052] d Dobs =d Dcal +Loss down

[0053] Get the downscaled observation record d Dobs Afterwards, the PIRNN-based inversion process can be performed on the downscaled model to update the parameters of the downscaled velocity model. The loss function used in the present invention is the mean square error function MSE, and the Adam optimizer is used for parameter update. The loss on the downscaled model is:

[0054] Loss=MSE(d Dobs ,d Dcal )

[0055] The calculated loss is back-propagated on the downscaled model to update its velocity parameter V. Finally, the updated velocity parameter V is converted to the original scale initial model through an upscaling operation. The above operation is the process of the next iterative inversion of the Dual-PIRNN network architecture. In this process, since the inversion operation of the present invention is performed on the downscaled model, the gradients that need to be saved in the back-propagation are fewer than those at the original scale, thereby reducing the computational cost.

[0056] 3. Scale decomposition

[0057] The purpose of scale conversion is to enable the information obtained by forward modeling on the original scale model to be used in the downscale model. This embodiment uses the Marmousi velocity model as an example and divides it into two different scales, the time scale down model and the space-time scale down model. The parameter setting example is shown in Table 1. Figure 5 As shown, the two scale models are called OriginalScale model and DownScale model. In this example, the grid number of the original scale model is 100x100, the grid spacing dx=dz=10, and the area size of 1000m x 1000m is simulated. The source function used is the sampling rate f s =4000Hz or f s =8000Hz Ricker wavelet, time point NT = 1600, the grid spacing dx, dz of the downscaled model is half of the original scale model, so the time point NT required for forward modeling is also half of the original scale, and the sampling rate of the source function used is f s =2000Hz, and the rest of the parameter settings are shown in Table 1.

[0058] Table 1. Marmousi model parameter settings at different scales

[0059]

[0060] It should be noted that in the above embodiment, the original-scale source sampling rate corresponding to time scaling is 4000 Hz, and the original-scale source sampling rate corresponding to spatiotemporal scaling is 8000 Hz. In actual applications, the parameter settings for time scaling and spatiotemporal scaling are determined based on the required reduction in resource consumption.

[0061] The technical effects of the present invention are described below with reference to specific data:

[0062] The Marmousi true velocity model and initial velocity model are established. The true velocity model is obtained by resampling the Marmousi public dataset, such as Figure 4 As shown in (a), the velocity range is [1500, 5500]. The initial velocity model is obtained by smoothing the real Marmousi model with different variance Gaussian. Figure 4 (b) and Figure 4 (c) The initial model settings when the variance σ = 5 and σ = 10 are respectively. The synthetic data set is calculated by forward propagation using the VSP observation system. Two downscaling strategies are used for inversion, namely time scale reduction and space-time scale reduction. The parameter settings are shown in Table 1. Adam is used as the optimization algorithm. The inversion results are shown in Figure 5 As shown by Figure 5It can be seen that the inversion results under the time scale can accurately identify the shallow velocity interface and also have a high resolution in the deep area. Figure 5 As shown in (a) and (c), good inversion results are obtained for initial models with different smoothness levels, demonstrating the robustness of the method. The inversion results of spatiotemporal downscaling can accurately predict the velocity boundaries in shallow areas and the velocity structure and trend in deep areas. Overall, due to the loss of more information in spatiotemporal decomposition, the effect is not as good as that of time downscaling, but it saves more GPU resources.

[0063] The speed curves at distances of 300m and 700m were selected as slices for comparison. The comparison results are as follows: Figure 6 As shown:

[0064] Figure 6 The three curves are given in the figure: the true value (True), the result value under the initial model of σ=10, and the result value under the initial model of σ=5. Figure 6 (a) and (b) show that the velocity values ​​in the shallow area at x = 300m and x = 700m in the high wrinkle area are very close to the true value, and at x = 700m the two are almost consistent. Under the inversion of the reduced time and space scale, Figure 6 As shown in (c) and (d), the velocity values ​​and changing trends can be accurately obtained in the shallow areas at x = 300 m and x = 700 m. There are periodic jumps in the deep areas, but the overall velocity trend is correct. The conclusion is consistent with the previous discussion. Figure 6 DownTime corresponds to down time, and DownTimeSpace corresponds to down time and space.

[0065] Table 2 shows the GPU resource usage under different inversion strategies. In the PIRNN full-waveform inversion architecture, the majority of GPU memory usage comes from storing the intermediate wavefield data related to the velocity parameter v. Using the Dual-PIRNN architecture, we reduce the number of numerical simulation time points and grid points using different downscaling strategies. Table 2 shows that the reduced-time inversion strategy reduces the intermediate activation usage from 32 GB to 16 GB, a reduction of approximately 50%. With the reduced-space-time inversion strategy, the intermediate activation usage is reduced from 32 GB to approximately 2 GB.

[0066] Table 2 GPU resource usage under different inversion strategies

[0067]

[0068] Model tests and real-world data demonstrate that the proposed Dual-PIRNN full-waveform inversion network structure significantly reduces the GPU resources required during the inversion process while maintaining a certain level of inversion accuracy. This demonstrates that the method effectively combines physical constraints with neural networks and scale decomposition, reducing computational costs while ensuring full-waveform inversion accuracy.

[0069] Those skilled in the art will appreciate that the embodiments described herein are intended to aid the reader in understanding the principles of the present invention, and it should be understood that the scope of the present invention is not limited to such specific descriptions and embodiments. Various modifications and variations are readily apparent to those skilled in the art. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention are intended to be included within the scope of the claims.

Claims

1. A dual-branch physical driven recurrent neural network VSP acoustic full waveform inversion method, characterized by: A dual-branch physics-driven recurrent neural network is constructed. The dual-branch physics-driven recurrent neural network specifically includes: an original-scale true velocity model forward model branch, an original-scale initial velocity model forward model branch, and a downscaled velocity model inversion branch. The original-scale true velocity model forward model branch, the original-scale initial velocity model forward model branch, and the downscaled velocity model inversion branch are each implemented using a physics-driven RNN forward network. The acoustic wave equation is solved by finite differences, and the velocity parameters in the equation are used as trainable parameters in the physics-driven RNN forward modeling network. The output of the forward modeling branch of the original-scale true velocity model is used as the observation value, and the output of the forward modeling branch of the original-scale initial velocity model is used as the prediction value. The loss between the prediction value and the observation value is input into the inversion branch of the downscaled velocity model to update the velocity parameters. After each iterative inversion, the velocity parameters on the downscaled velocity model inversion branch are updated by upsampling to the forward model branch of the original-scale initial velocity model for another forward modeling, and the update cycle is repeated. The loss between the predicted value and the observed value is recorded as Loss original ; The average pooling-based downsampling method is used to Loss original Extract and get the loss that matches the downscaled velocity model Loss down ; In the downscaled velocity model inversion branch, the VSP shot gather records on the downscaled velocity model are obtained by forward modeling based on the physics-driven RNN forward modeling network. d Dcal ; Will Loss original and d Dcal Add together and get the observation record on the downscaled velocity model d Dobs .

2. The VSP acoustic full waveform inversion method of a dual-branch physical driven recurrent neural network according to claim 1 is characterized in that: The input of the physics-driven RNN forward modeling network is the velocity parameters, source location, and wave field information at the previous moment, and the output is the VSP shot gather record at the current moment.

3. The VSP acoustic full waveform inversion method of a dual-branch physical driven recurrent neural network according to claim 2 is characterized in that: Loss original The loss is obtained by calculating the L1 norm between the predicted value and the observed value.

4. The VSP acoustic full waveform inversion method of a dual-branch physical driven recurrent neural network according to claim 3 is characterized in that: The loss function used in the downscaling velocity model is the mean square error function.

5. The VSP acoustic full waveform inversion method of a dual-branch physical driven recurrent neural network according to claim 4 is characterized in that: The downscaling velocity model is specifically a downscaling time scale model or a downscaling space-time scale model.

Images