Seismic Velocity Inversion via Physical Embedding Recurrent Neural Networks (RNN)

Abstract

Seismic velocity inversion is one of the most critical issues in the field of seismic exploration and has long been the focus of numerous experts and scholars. In recent years, the advancement of machine learning technologies has infused new vitality into the research of seismic velocity inversion and yielded a wealth of research outcomes. Typically, seismic velocity inversion based on machine learning lacks control over physical processes and interpretability. Starting from wave theory and the physical processes of seismic data acquisition, this paper proposes a method for seismic velocity model inversion based on Physical Embedding Recurrent Neural Networks. Firstly, the wave equation is a mathematical representation of the physical process of acoustic waves propagating through a medium, and the finite difference method is an effective approach to solving the wave equation. With this in mind, we introduce the architecture of recurrent neural networks to describe the finite difference solution of the wave equation, realizing the embedding of physical processes into machine learning. Secondly, in seismic data acquisition, the propagation of acoustic waves from multiple sources through the medium represents a high-dimensional causal time series (wavefield snapshots), where the influential variable is the velocity model, and the received signals are the observations of the wavefield. This forms a forward modeling process as the forward simulation of the wavefield equation, and the use of error back-propagation between observations and calculations as the velocity inversion process. Through time-lapse inversion and by incorporating the causal information of wavefield propagation, the non-uniqueness issue in velocity inversion is mitigated. Through mathematical derivations and theoretical model analyses, the effectiveness and rationality of the method are demonstrated. In conjunction with simulation results for complex models, the method proposed in this paper can achieve velocity inversion in complex geological structures.

1. Introduction

Seismic waveform inversion is an essential geophysical method that enables the inversion of physical parameters of the subsurface medium, such as velocity, density, and elastic moduli, from seismic data [1]. These parameters are of great significance in fields such as oil and gas exploration, earthquake hazard assessment, and geological structure analysis [2]. However, seismic waveform inversion is also a non-linear, multi-modal, non-convex optimization problem and faces various challenges, such as dependency on initial models, susceptibility to local minima traps, and high computational complexity [3].

To overcome these challenges, deep learning techniques have been introduced to seismic inversion. Röth and Tarantola were the first to propose the use of neural networks to transform time-domain seismic data into acoustic velocity models in 1994 [4]. Nath et al., in 1999, trained neural networks with synthetic data, resulting in networks capable of estimating formation velocities using inter-well data [5,6]. In 2019, Yang and Ma developed a fully convolutional neural network for inverting P-wave velocity models from raw seismic data [7]. A commonality among these networks is the need for extensive labeled data for training, and once trained, the neural networks can be used as inversion solvers. However, in the field of seismic exploration, acquiring a large labeled dataset is challenging.

Physics-Embedded Neural Networks (PINNs) offer a novel approach to seismic velocity modeling [8]. The fundamental idea behind PINNs is to use neural networks to approximate unknown functions or parameters and incorporate the residuals of physical equations as part of the loss function. The neural networks are trained by minimizing this loss function, resulting in solutions that satisfy physical constraints [9]. Moseley et al. (2020) utilized PINNs as solvers for acoustic wave forward propagation [10,11], while Smith et al. [12] and Waheed et al. [13] applied PINNs to the Eikonal equation for first-arrival prediction and travel-time tomography, respectively. Song et al. employed PINNs to solve the frequency-domain anisotropic acoustic wave equation [14]. Additionally, PINNs have been successfully used for forward and inverse analyses of various complex systems in cases with limited sensor measurements [15,16,17,18,19,20] or without labeled data [21,22,23,24,25].

The recurrent neural network, as a special kind of neural network, introduces a time dependence, which is more compatible with the wave field propagation process. In this paper, a recurrent neural network is used to simulate the two-dimensional time-domain finite-difference acoustic wave equation as the forward process using a similar idea of full-waveform inverse forward simulation combined with backpropagation correction, and the velocity is used as a network parameter [26]. The velocity model is updated during the backward propagation of the network. This approach enables the integration of information from the acoustic wave equation into the network, guiding the optimization direction of the network, and avoiding explicit solutions for wavefields and gradients, thereby reducing computational costs. Furthermore, it leverages the non-linear fitting capabilities of neural networks to enhance inversion accuracy and robustness. Additionally, given the temporal causality in wavefields, where the current computed wavefield is related to previous wavefields, this paper adopts a time-lapse inversion approach. This ensures that the residuals of the wavefields are minimized at each timestep, guaranteeing the accuracy of subsequent wavefield computations and alleviating the non-uniqueness issue in velocity inversion. The effectiveness and superiority of the method proposed in this paper are verified through numerical experiments.

2. Materials and Methods

2.1. Finite-Difference Acoustic Wave Equation

The numerical simulation of wave equation forward modeling is fundamental to inversion. The time-space domain finite difference algorithm is widely used in geophysical forward modeling due to its simplicity and high computational efficiency, and it is employed to simulate the propagation of seismic waves. In the process of solving partial differential wave equations using the finite difference method, the computation area is first discretized into a grid, then grid differences are utilized to approximate the derivative operators, and the wave equation is converted into a difference equation.

In a 2D isotropic acoustic medium, the time-domain wave equation is

(1)${\nabla }^2\mathrm{u}\left(\mathrm{r},\mathrm{t}\right)=\frac{1}{{\mathrm{v}}^2\left(\mathrm{r}\right)}\cdot \frac{\partial {\mathrm{u}}^2\left(\mathrm{r},\mathrm{t}\right)}{\partial {\mathrm{t}}^2}+\mathrm{s}\left(\mathrm{r},\mathrm{t}\right)\mathsf{\delta }\left(\mathrm{r}-{\mathrm{r}}_{\mathrm{s}}\right)$

where $\nabla {\mathrm{u}}^2\left(\mathrm{r},\mathrm{t}\right)$ represents the pressure or displacement at coordinate r at time t, ${\mathrm{V}}^2\left(\mathrm{r}\right)$represents the speed of sound at coordinate r, $\mathrm{s}\left(\mathrm{r},\mathrm{t}\right)$ represents the source function at coordinate r at time t, and $\mathsf{\delta }\left(\mathrm{r}-{\mathrm{r}}_{\mathrm{s}}\right)$ is a Dirac (Delta) function, indicating that energy is released only at the source position.

For the 2D medium, spatial discretization is undertaken by converting coordinates (x, y) into grid points (i, j), where i represents the index in the x direction and j represents the index in the y direction. Uniform grid division is usually adopted, dividing the space into a series of equal-sized cells. Time t is discretized into time steps Δt, allowing time to progress in fixed intervals. Based on discrete spatial grids and time steps, difference approximations are used to replace derivative terms. A centered difference approximation is used to approximate the second-order spatial derivative $\nabla {\mathrm{u}}^2\left(\mathrm{r},\mathrm{t}\right)$:

(2)$\nabla {\mathrm{u}}^2\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)\approx \frac{\mathrm{u}\left(\mathrm{i}+1,\mathrm{j},\mathrm{t}\right)-2\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)+\mathrm{u}\left(\mathrm{i}-1,\mathrm{j},\mathrm{t}\right)}{\mathsf{\Delta }{\mathrm{x}}^2}+\frac{\mathrm{u}\left(\mathrm{i},\mathrm{j}+1,\mathrm{t}\right)-2\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)+\mathrm{u}\left(\mathrm{i},\mathrm{j}-1,\mathrm{t}\right)}{\mathsf{\Delta }{\mathrm{y}}^2}$

where Δx and Δy represent the grid spacing in the x and y directions, respectively.

A second-order centered difference approximation is used to approximate the second-order time derivative $\partial {\mathrm{u}}^2\left(\mathrm{r},\mathrm{t}\right)/\partial {\mathrm{t}}^2$:

(3)$\frac{\partial {\mathrm{u}}^2\left(\mathrm{r},\mathrm{t}\right)}{\partial {\mathrm{t}}^2}\approx \frac{\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}+\mathsf{\Delta }\mathrm{t}\right)-2\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)+\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}-\mathsf{\Delta }\mathrm{t}\right)}{\mathsf{\Delta }{\mathrm{t}}^2}$

Substituting the difference approximations into the original wave equation yields the 2D time-domain finite difference acoustic wave equation:

(4)$\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}+\mathsf{\Delta }\mathrm{t}\right)={\mathrm{V}}^2\left(\mathrm{i},\mathrm{j}\right)\mathsf{\Delta }{\mathrm{t}}^2\left[\nabla {\mathrm{u}}^2\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)-\mathrm{s}\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)\mathsf{\delta }\left(\mathrm{i},\mathrm{j}-{\mathrm{r}}_{\mathrm{s}}\right)\right]+2\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}\right)-\mathrm{u}\left(\mathrm{i},\mathrm{j},\mathrm{t}-\mathsf{\Delta }\mathrm{t}\right)$

Through the difference acoustic wave equation, the wavefield at the current time can be calculated using the wavefields of the past two time steps, the source, and the velocity model. This relationship can be seen as a function g:

(5)$\mathrm{u}\left(\mathrm{r},\mathrm{t}+\mathsf{\Delta }\mathrm{t}\right)=\mathrm{g}\left(\mathrm{u}\left(\mathrm{r},\mathrm{t}\right),\mathrm{u}\left(\mathrm{r},\mathrm{t}-\mathsf{\Delta }\mathrm{t}\right),\mathrm{s}\left(\mathrm{r},\mathrm{t}\right),{\mathrm{v}}^2\left(\mathrm{r}\right)\right)$

The wavefield at each moment can be considered as part of the time series data. The calculation process of the difference acoustic wave equation is shown in Figure 1 and Figure 2. A more in-depth study of the propagation process of the acoustic wave equation in more complex media was carried out by Abolghasmifar et al. [27].

2.2. Recurrent Neural Network with Acoustic Wave Equation Embedding

The structure of the RNN (Recurrent Neural Network) is shown in Figure 3. By comparison, it can be observed that the propagation structure of the difference wave equation has high similarity with the network structure of the RNN, and existing experiments indicate that the dynamics of the wave equation are conceptually equivalent to the RNN. Therefore, the acoustic difference equation can be embedded into the RNN structure, which enables a better understanding of the propagation characteristics of waves and maps them into the neural network structure [28].

Considering the dependence of the acoustic difference equation on the wavefields of the two previous time steps, these wavefields can be used as hidden states in the RNN. The source function at the current time step can be used as the input, and the signal at the detector position for each time step ‘t’ can be used as the shot gather at time ‘t’. The velocity model in the difference equation can be treated as a trainable parameter within the network. By using the error between the observed and predicted values for backward correction, this process can be viewed as the inversion of the velocity model.

This approach allows the information from the acoustic wave equation to be incorporated into the network, guiding the optimization direction of the network and avoiding the explicit calculation of gradients, which reduces computational costs. Additionally, the non-linear fitting capabilities of neural networks can be leveraged to improve inversion accuracy and robustness.

The structure of the RNN with embedded difference acoustic wave equation and the structure of operators are shown in Figure 4 and Figure 5, respectively.

2.3. Time-by-Time Inversion Algorithm

In RNNs, two common network updating algorithms are generally employed: (1) Compute the loss and update the network at each time step: At each time step, the loss is computed and used to update the network parameters. This is a traditional time-based loss update method. (2) Compute the loss at each time step, sum them, and then update the network: At each time step, the current loss is computed and accumulated with losses from previous time steps. The network parameters are then updated using the accumulated sum of losses. This method integrates information from previous time steps into the current time step loss to provide global information.

In the RNN embedded with the acoustic wave equation, updating the network by computing the loss at each time step would result in changing the velocity model at different time steps. According to the calculation formula of the difference wave Equation (5), the velocity needs to remain consistent throughout the process of solving the acoustic wave equation. Therefore, computing the loss and updating the parameters separately at each time step would violate the constraints of the acoustic wave equation, leading to the inversion not converging. This article introduces a new network loss update method called the step-by-step inversion algorithm. Considering the causality of wavefield propagation, the calculation of the wavefield at the current time step is related to the previous wavefields. If errors are present in the calculation of the preceding wavefields, this will introduce additional errors in the subsequent wavefield calculations. Therefore, it is necessary to ensure the correctness of the preceding wavefield calculations before calculating subsequent wavefields. The step-by-step inversion adopts the idea of ensuring the correctness of the preceding wavefields for inversion. Specifically, starting from the beginning of wavefield propagation, the sum of the losses of the wavefield at each time step and its preceding time steps is computed to calculate the total loss. The model parameters are updated according to the total loss to ensure the correctness of the wavefield at the current time step. Compared to the traditional method of independently calculating the loss and updating the network at each time step or calculating the sum of the losses for all time steps, this method has significant advantages. Firstly, step-by-step inversion not only minimizes the overall residual but also ensures that the residual of the wavefield at each time step is as small as possible. In addition, by ensuring that the residual of the wavefield at the previous time step is small before propagating to the next time step, this method incorporates causal information of wavefield propagation, satisfying the propagation conditions of the acoustic wave equation. This avoids errors at the current time step that could be caused by errors from the previous time step, improving the accuracy of the inversion. By using a velocity model that has been partially corrected through inversion as the network parameters for the next time step, the initial model for the next time step inversion is altered, which can alleviate the dependency on the initial model in velocity inversion.

The mean squared error loss is calculated at each time step by comparing the output with the label:

(6)$\mathsf{\delta }{\mathrm{d}}_{\mathrm{t}}={{(\mathrm{d}}_{\mathrm{t}}-{{\mathrm{d}}^{\prime }}_{\mathrm{t}})}^2$

where ${\mathrm{d}}_{\mathrm{t}}$ is the output of the network at time t, and ${{\mathrm{d}}^{\prime }}_{\mathrm{t}}$ is the label sequence at time t. The residual at time t used for computing the gradient through backpropagation is:

(7)$\mathrm{l}\mathrm{o}\mathrm{s}{\mathrm{s}}_{\mathrm{t}}=\sum _{{\mathrm{r}}_{\mathrm{s}}}\sum _0^{\mathrm{t}}\mathsf{\delta }{\mathrm{d}}_{\mathrm{t}}$

where ${\mathrm{r}}_{\mathrm{s}}$ represents different shot gathers. The corrections to the model parameters can be derived from the derivatives of the total losses to the parameters:

(8)$\mathrm{g}\left(\mathrm{r}\right)=\frac{\partial \mathrm{l}\mathrm{o}\mathrm{s}{\mathrm{s}}_{\mathrm{t}}}{\partial \mathrm{V}\left(\mathrm{r},\mathrm{t}\right)}$

(9)$\mathrm{g}\left(\mathrm{r}\right)=\sum _0^{\mathrm{t}}\frac{\mathrm{l}\mathrm{o}\mathrm{s}{\mathrm{s}}_{\mathrm{t}}}{\partial \mathrm{u}\left(\mathrm{r},\mathrm{t}\right)}\cdot \frac{\partial \mathrm{u}\left(\mathrm{r},\mathrm{t}\right)}{\partial \mathrm{V}\left(\mathrm{r}\right)}$

where the partial derivatives of the loss with respect to the wave field can be computed by Equation (4) and the backpropagation of the network is accomplished by passing the computation to the initial moment through the chain rule. After obtaining the parameter update quantity, the parameters are updated by Equation (10), which show the learning rate.

(10)${\mathrm{V}}_{\mathrm{n}+1}\left(\mathrm{r}\right)={\mathrm{V}}_{\mathrm{n}}\left(\mathrm{r}\right)-ɑ\cdot \mathrm{g}\left(\mathrm{r}\right)$

The gradient is propagated through the chain rule, so the velocity correction is related to ${\mathrm{d}}_{\mathrm{t}}$, which causes the residual. Through the two-dimensional finite difference acoustic wave equation in the time domain (5), it is known that ${\mathrm{d}}_{\mathrm{t}}$ is related to the velocity model and the wavefields of the previous two time steps. The spatial domain variations at different time steps are caused by the difference in the wavefields of the previous time step. Therefore, the range of velocity correction gradually expands with the range of the wavefield. Through the step-by-step correction algorithm in the time domain, the velocity model can be gradually corrected in the spatial domain, starting from the source and progressively correcting the velocity model in all directions. The gradual correction in the spatial domain makes the inversion parameters smaller, thereby reducing the non-uniqueness of the inversion problem. Moreover, a more accurate shallow velocity can improve the accuracy and stability of deep velocity corrections. The schematic diagram of the step-by-step inversion process is shown in Figure 6.

3. Results

To validate the effectiveness of the Step-by-Step Inversion Algorithm and its advantages over the Sum-of-Losses Update method, numerical experiments were conducted using both a Homogeneous Layer Model and a Sloping Layer Model for comparison and analysis. Ultimately, the Step-by-Step Inversion Algorithm was applied to synthetic data computed using the Marmousi model to verify the feasibility of this algorithm.

3.1. Homogeneous Layer Model

In order to validate the efficacy of the Step-by-Step Inversion Algorithm, a Homogeneous Layer Model was constructed. In this model, the velocities from the shallow to deep regions are set as [3, 3.5, 4] km/s. Synthetic datasets were computed through forward propagation using a Vertical Seismic Profile (VSP) observation system, Orthogonal parameter points are 100 × 100, with a spatial step of 10 m and a time step of 0.5 ms. The initial model has velocities set as [3, 3.4, 3.9] km/s from shallow to deep regions, simulating a scenario where the initial model is relatively accurate. The inversion was performed using the Adam optimization algorithm, with identical parameter settings for all tests. The Homogeneous Layer Model and its inversion results are depicted in Figure 7.

As can be seen from Figure 7, the Step-by-Step Inversion Algorithm yields smoother results both in shallow and deep areas, especially near the source points close to the surface and the lateral regions. Slices were taken at depths of 500 m and 800 m for velocity profile comparisons, with the results shown in Figure 8 and Figure 9, respectively.

It can be deduced from Figure 8 and Figure 9 that, in terms of inversion accuracy, the Step-by-Step Inversion generally achieves higher accuracy compared to the Sum-of-Losses Inversion. The inversion results are closer to the true model. Moreover, by employing a time-stepped inversion, which takes into account the propagation of waves, the accuracy of the preceding wavefield is ensured. This, in turn, leads to more precise calculations in subsequent wavefield propagations, manifesting in the inversion results as more accurate representations on the lateral sides and deeper regions of the model.

3.2. Sloping Layer Model

To evaluate the performance of the Step-by-Step Inversion Algorithm in inclined geological structures, as well as its weak dependence on the initial model, a sloping layer model was constructed. In this second model, the velocities increase from the shallow to deeper layers, following [2.5, 3.5, 4.5] km/s, respectively. The synthetic dataset was computed using a VSP (Vertical Seismic Profile) observation system via forward modeling, and the orthogonal parameters are consistent with model I. The initial model had velocities [2, 3, 4] km/s, with significant discrepancies between the initial and true models. This was to validate the weak dependency of the Step-by-Step Inversion Algorithm on the initial model, demonstrating that accurate inversion is still achievable despite large errors in the initial model. The inversion employed the Adam optimization algorithm with identical parameter settings. The results of the second model and its inversion are shown in Figure 10.

As can be observed in Figure 10, despite the substantial error in the initial model in the context of a sloping layer model, the Step-by-Step Inversion Algorithm exhibited superior performance overall. This suggests that the algorithm has a lower dependency on the initial model and is capable of achieving reasonably accurate inversions even when the initial model has large discrepancies. In contrast, the traditional Sum-of-Losses Inversion performed poorly, indicating a higher dependency on the initial model. Velocity profiles at depths of 150 m, 500 m, and 800 m were selected for comparison and are displayed in Figure 11, Figure 12 and Figure 13, respectively.

Based on the cross-sectional velocity comparison plots, it can be observed that, even with significant errors in the initial model, the Step-by-Step Inversion Algorithm was able to more accurately invert velocities at various depths, including 150 m, 500 m, and 800 m, demonstrating superior precision in the results.

3.3. Marmousi Model

To verify the effectiveness of the algorithm under more complex real-world conditions, the experiments used a Marmousi model. The model was created through a collaboration of numerous geologists, geophysicists, and computer scientists and is capable of simulating actual complex geological conditions and is geologically consistent [29]. A relatively accurate and smooth velocity model was used as the initial model for the inversion. The real and initial models are shown in Figure 14. The VSP observation system was used to generate synthetic data from the real model with parameters of 100 × 200, a spatial step of 10 m, and a time step of 0.5 ms.

The Marmousi model and inversion results are shown in Figure 15. As can be seen from the inversion results, both the overall inversion and time-by-time inversion results are relatively close to the real model. A comparison of horizontal velocity values at depths of 500 m and 850 m is shown in Figure 16 and Figure 17, and a comparison of depth velocity values at a horizontal distance of 1000 m is shown in Figure 18. From the velocity comparison curve, it can be seen that both the time-by-time inversion and the overall inversion are relatively close to the real model. Moreover, the velocity of the time-by-time inversion results is smoother, with fewer velocity anomalies, making the inversion more stable.

4. Discussion

We believe that maximizing the use of known information in the velocity inversion process can improve the accuracy of the inversion. The causality of wavefield propagation is introduced in physically embedded recurrent neural networks. The recurrent neural network can characterize the forward propagation of the wave field well, and the time-by-time inversion can strengthen the causality of the wave field propagation. In time-by-time inversion, the inversion result of the previous moment is used as the initial model of the inversion in the subsequent moment, which can make the wave field propagation not only causal but also correct. Through the experiments and comparisons of the three different models, it can be seen that when the initial model is better, the time-by-time inversion can achieve more accurate results than the traditional overall inversion, which shows that it is effective to ensure the accuracy of causality. When the initial model is poor, the time-by-time inversion is still more effective, which indicates that the time-by-time inversion is less dependent on the initial model. Finally, in the Marmousi model, the time-by-time inversion achieved more accurate results on both sides and in the deeper regions. It is worth mentioning that since the time-by-time inversion has to ensure the correctness of the cause and effect, when the error of the initial model is large, the time-by-time inversion needs more time to correct the wavefield at each moment. The study by Rointan et al. provides guidance on initial model selection [30].

5. Conclusions

This paper presents a physically embedded recurrent neural network (RNN)-based time-by-time inversion algorithm. Using the time-by-time inversion algorithm, the causal information in the time evolution of the wave equation is mapped into spatial information for velocity model correction. The velocity model is corrected in the spatial region according to the wave field propagation range. The corrected model will be used as the initial model for subsequent time steps covering a larger time range, thus reducing the dependence of the inversion on the initial model. In addition, controlling the correction region can reduce the inversion parameters and effectively reduce the non-uniqueness of the inversion solution within a single time step, and eventually reduce the overall non-uniqueness, while improving the accuracy and stability of the inversion. Overall, the time-by-time inversion can achieve more accurate results than the overall inversion, but it may require longer computation time in the case of large initial model errors. The effectiveness and superiority of the method proposed in this paper are verified by numerical experiments.

Author Contributions

Conceptualization, C.L. and C.Z.; methodology, C.L.; software, C.Z.; validation, C.Z. and C.L.; formal analysis, C.L.; investigation, C.L.; resources, C.Z.; data curation, C.Z.; writing—original draft preparation, C.Z.; writing—review and editing, C.L.; visualization, C.L.; supervision, C.L.; project administration, C.L.; funding acquisition, C.L. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data are contained within the article.

Acknowledgments

We are grateful to the reviewers for their valuable comments.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures

Figure 1. Schematic of wave field propagation.

Figure 2. Schematic diagram of differential acoustic wave equation.

Figure 3. Schematic diagram of RNN network.

Figure 4. Schematic diagram of the structure of the RNN network embedded in the acoustic wave equation.

Figure 5. Embedding of the acoustic wave equation into the RNN operator.

Figure 6. Schematic diagram of the time-by-time inversion algorithm.

Figure 7. Homogeneous layer model and inversion results. (a) True model. (b) Total loss inversion. (c) Time by time inversion.

Figure 8. Comparison of 500 m deep velocity profiles of Homogeneous Layer Model.

Figure 9. Comparison of 800 m deep velocity profiles of Homogeneous Layer Model.

Figure 10. Sloping layer model and inversion results. (a) True model. (b) Total loss inversion. (c) Time by time inversion.

Figure 11. Comparison of 150 m deep velocity profiles of Sloping Layer Model.

Figure 12. Comparison of 500 m deep velocity profiles of Sloping Layer Model.

Figure 13. Comparison of 800 m deep velocity profiles of Sloping Layer Model.

Figure 14. Marmousi model and initial model. (a) True model. (b) Initial model.

Figure 15. Marmousi model and inversion results. (a) True model. (b) Total loss inversion. (c) Time by time inversion.

Figure 16. Comparison of 500 m deep velocity profiles of Marmousi Model.

Figure 17. Comparison of 850 m deep velocity profiles of Marmousi Model.

Figure 18. Comparison of 1000 m distance velocity profiles of Marmousi Model.

References

1. Liu, Y.; He, B.; Zheng, Y. Controlled-order multiple waveform inversion. Geophysics; 2020; 85, pp. R243-R250. [DOI: https://dx.doi.org/10.1190/geo2019-0658.1]

2. Kuang, W.; Yuan, C.; Zhang, J. Real-time determination of earthquake focal mechanism via deep learning. Nat. Commun.; 2021; 12, 1432. [DOI: https://dx.doi.org/10.1038/s41467-021-21670-x]

3. Zhang, J.; Zhang, H.; Chen, E.; Zheng, Y.; Kuang, W.; Zhang, X. Real-time earthquake monitoring using a search engine method. Nat. Commun.; 2014; 5, 5664. [DOI: https://dx.doi.org/10.1038/ncomms6664] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25472861]

4. Roth, G.; Tarantola, A. Neural Networks and Inversion of Seismic Data. J. Geophys. Res. Solid Earth; 1994; 99, pp. 6753-6768. [DOI: https://dx.doi.org/10.1029/93JB01563]

5. Nath, S.K.; Vasu, R.; Pandit, M. Wavelet based compression and denoising of optical tomography data. Opt. Commun.; 1999; 167, pp. 37-46. [DOI: https://dx.doi.org/10.1016/S0030-4018(99)00294-1]

6. Kumar, N.S.; Subrata, C.; Kumar, S.S.; Nilanjan, G. Velocity inversion in cross-hole seismic tomography by counter-propagation neural network, genetic algorithm and evolutionary programming techniques. Geophys. J. Int.; 1999; 138, pp. 108-124.

7. Yang, F.; Ma, J. Deep-learning inversion: A next-generation seismic velocity model building method. Geophysics; 2019; 84, pp. R583-R599. [DOI: https://dx.doi.org/10.1190/geo2018-0249.1]

8. Ren, P.; Rao, C.; Sun, H.; Liu, Y. Physics-informed neural network for seismic wave inversion in layered semi-infinite domain. Comput. Methods Appl. Math.; 2022; 22, pp. 1-24.

9. Raissi, M.; Perdikaris, P.; Karniadakis, G.E. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations. J. Comput. Phys.; 2019; 378, pp. 686-707. [DOI: https://dx.doi.org/10.1016/j.jcp.2018.10.045]

10. Moseley, B.; Nissen-Meyer, T.; Markham, A. Deep learning for fast simulation of seismic waves in complex media. Solid Earth; 2020; 11, pp. 1527-1549. [DOI: https://dx.doi.org/10.5194/se-11-1527-2020]

11. Moseley, B.; Markham, A.; Nissen-Meyer, T. Solving the wave equation with physics-informed deep learning. arXiv; 2020; arXiv: 2006.11894

12. Smith, J.D.; Azizzadenesheli, K.; Ross, Z.E. Eikonet: Solving the eikonal equation with deep neural networks. IEEE Trans. Geosci. Remote Sens.; 2020; 59, pp. 10685-10696. [DOI: https://dx.doi.org/10.1109/TGRS.2020.3039165]

13. Waheed, U.B.; Alkhalifah, T.; Haghighat, E.; Song, C.; Virieux, J. PINNtomo: Seismic tomography using physics-informed neural networks. arXiv; 2021; arXiv: 2104.01588

14. Song, C.; Alkhalifah, T.; Waheed, U.B. Solving the frequency-domain acoustic VTI wave equation using physics-informed neural networks. Geophys. J. Int.; 2021; 225, pp. 846-859. [DOI: https://dx.doi.org/10.1093/gji/ggab010]

15. Karniadakis, G.E.; Kevrekidis, I.G.; Lu, L.; Perdikaris, P.; Wang, S.; Yang, L. Physics-informed machine learning. Nat. Rev. Phys.; 2021; 3, pp. 422-440. [DOI: https://dx.doi.org/10.1038/s42254-021-00314-5]

16. Rao, C.; Sun, H.; Liu, Y. Physics-informed deep learning for incompressible laminar flows. Theor. Appl. Mech. Lett.; 2020; 10, pp. 207-212. [DOI: https://dx.doi.org/10.1016/j.taml.2020.01.039]

17. Raissi, M.; Yazdani, A.; Karniadakis, G.E. Hidden fluid mechanics: Learning velocity and pressure fields from flow visualizations. Science; 2020; 367, pp. 1026-1030. [DOI: https://dx.doi.org/10.1126/science.aaw4741]

18. Lu, L.; Meng, X.; Mao, Z.; Karniadakis, G.E. Deepxde: A deep learning library for solving differential equations. SIAM Rev.; 2021; 63, pp. 208-228. [DOI: https://dx.doi.org/10.1137/19M1274067]

19. Haghighat, E.; Raissi, M.; Moure, A.; Gomez, H.; Juanes, R. A physics-informed deep learning framework for inversion and surrogate modeling in solid mechanics. Comput. Methods Appl. Mech. Eng.; 2021; 379, 1137411. [DOI: https://dx.doi.org/10.1016/j.cma.2021.113741]

20. Ren, P.; Rao, C.; Liu, Y.; Ma, Z.; Wang, Q.; Wang, J.-X.; Sun, H. Physics-informed deep super-resolution for spatiotemporal data. arXiv; 2022; arXiv: 2208.01462[DOI: https://dx.doi.org/10.1016/j.jcp.2023.112438]

21. Sun, L.; Gao, H.; Pan, S.; Wang, J.-X. Surrogate modeling for fluid flows based on physics-constrained deep learning without simulation data. Comput. Methods Appl. Mech. Eng.; 2019; 361, 112732. [DOI: https://dx.doi.org/10.1016/j.cma.2019.112732]

22. Rao, C.; Sun, H.; Liu, Y. Physics-informed deep learning for computational elastodynamics without labeled data. J. Eng. Mech.; 2021; 147, 04021043. [DOI: https://dx.doi.org/10.1061/(ASCE)EM.1943-7889.0001947]

23. Gao, H.; Sun, L.; Wang, J.-X. Phygeonet: Physics-informed geometry-adaptive convolutional neural networks for solving parameterized steady-state pdes on irregular domain. J. Comput. Phys.; 2021; 428, 110079. [DOI: https://dx.doi.org/10.1016/j.jcp.2020.110079]

24. Ren, P.; Rao, C.; Liu, Y.; Wang, J.-X.; Sun, H. Phycrnet: Physics-informed convolutional-recurrent network for solving spatiotemporal pdes. Comput. Methods Appl. Mech. Eng.; 2022; 389, 114399. [DOI: https://dx.doi.org/10.1016/j.cma.2021.114399]

25. Gao, H.; Zahr, M.J.; Wang, J.-X. Physics-informed graph neural Galerkin networks: A unified framework for solving PDE-governed forward and inverse problems. Comput. Methods Appl. Mech. Eng.; 2022; 390, 114502. [DOI: https://dx.doi.org/10.1016/j.cma.2021.114502]

26. Alaei, N.; Soleimani Monfared, M.; Roshandel Kahoo, A.; Bohlen, T. Seismic Imaging of Complex Velocity Structures by 2D Pseudo-Viscoelastic Time-Domain Full-Waveform Inversion. Appl. Sci.; 2022; 12, 7741. [DOI: https://dx.doi.org/10.3390/app12157741]

27. Abolghasemifar, S.A.; Monfared, M.S.; Aghajani, H. Acoustic wave propagation in strongly heterogeneous media: A lattice boltzmann approach. Int. Multidiscip. Sci. GeoConference SGEM; 2009; 1, pp. 643-648.

28. Sun, J.; Niu, Z.; Innanen, K.A.; Li, J.; Trad, D.O. A theory-guided deep-learning formulation and optimization of seismic waveform inversion. Geophysics; 2020; 85, pp. R87-R99. [DOI: https://dx.doi.org/10.1190/geo2019-0138.1]

29. Brougois, A.; Bourget, M.; Lailly, P.; Poulet, M. Marmousi, model and data. EAEG Workshop-Practical Aspects of Seismic Data Inversion; European Association of Geoscientists & Engineers: Utrecht, The Netherlands, 1990; cp-108-00002

30. Rointan, A.; Soleimani Monfared, M.; Aghajani, H. Improvement of seismic velocity model by selective removal of irrelevant velocity variations. Acta Geod. Et Geophys.; 2021; 56, pp. 145-176. [DOI: https://dx.doi.org/10.1007/s40328-020-00329-x]