Nonlinear amplitude versus angle inversion using hybrid quantum ant colony optimization and the exact Zoeppritz equation

Prestack amplitude versus offset (AVO) inversion is an essential tool to estimate elastic properties. Approximations of P-wave reflection coefficients based on the Zoeppritz equation are limited by the assumption of weak contrast and are inaccurate for far offset, and the existing nonlinear inversion method, like ant colony optimization (ACO), cannot converge to the global optimal solution within a limited time. Consequently, we propose a nonlinear AVO inversion method based on the exact Zoeppritz equation using hybrid quantum ant colony optimization (HQACO). The Zoeppritz equation is a forward modeling approach with high accuracy at far offset and is valid at the strong-contrast interface. HQACO, as a global optimization algorithm, not only has the potential of ACO to search global solutions but also exploits the power of quantum computing to speed up optimization procedures. The self-adaptive rotating strategy is proposed to improve the flexibility and efficiency of the conventional quantum rotating gate largely. Moreover, the quantum gate is introduced to enhance the global search ability. Numerical results show that the average evaluation number required for the global solution decreases largely and the solution converges to the global optimal solution closely using the proposed method. Synthetic applications verify that HQACO shows the most reliable exploring competence to find the global solution. The new inversion method is not only suitable for the strong-contrast interface with noise but also shows promising accuracy using the wide offset range.

1 Introduction

Prestack amplitude variation with offset (AVO) inversion is one of the commonly used techniques for estimating reservoir elastic parameters (Ostrander, 1984). In AVO inversion, an objective function is usually designed to fit the prior constraint and observation data by finding the minimum value of the objective function to obtain the true solution of the problem. Moreover, seismic data are band limited and always perturbed by noise, which inevitably leads to non-uniqueness of the solution, that is, ill-posed problem. In this sense, the inversion problem is a trade-off optimization problem with non-unique solutions and high-dimensional nonlinear characteristics. In order to improve the computational efficiency and stabilize inversion, many approximations have been developed to eliminate the nonlinear relationship between the seismic response and the elastic parameters by linearization under the assumptions of weak contrast and near offset limitations (Aki and Richards, 2002; Shuey, 1985; Smith and Gidlow, 1987; Zong et al., 2013; Zhang and Li, 2013; Lu et al., 2015; Zhang et al., 2019; Cheng et al., 2022a). The global optimization algorithm, however, can never reach the true model parameters under such assumption without considering the computational efficiency (Lu et al., 2015). Consequently, an exact inversion equation, like the Zoeppritz equation, is a necessary condition for AVO inversion, and an accurate inversion method with acceptable computation efficiency facilitates the inverse problem.

Many inversion methods are proposed to solve the nonlinear inverse problem (Stoffa and Sen, 1991; Li and Mallick, 2015; Sen and Stoffa, 2013; Mallick and Adhikari, 2015; Aleardi and Mazzotti, 2017; Zhou et al., 2017; Aleardi et al., 2019; Liu et al., 2022). The linear inversion method, that is, deterministic algorithm, solves the nonlinear function between model parameters and observational data as a linear problem and iteratively modifies the initial model until the convergence conditions are met. This inversion method has high calculation efficiency, requires less memory, and is concise and commonly used in practice (Buland and Omre, 2003; Russell et al., 2011; Zhang et al., 2019). However, the inversion result depends heavily on the initial model, and it may converge to a local solution rather than a global solution, which may cause large deviation between the solution and the actual model (Liu et al., 2022; Liu et al., 2023). Considering that the exact equations are high dimensional and nonlinear, extra assumptions like linearization make global optimization methods miss the global optimal solution, which means the desired global optimization methods still cannot acquire the real solutions without an exact forward modeling algorithm. Therefore, the inversion method based on the exact Zoeppritz equation can meet the requirements of the global optimization methods. Although the most ideal inversion method in such a problem is the exhaustive search method, it is also the most unrealistic method because the optimization method requires infinite calculation time whatever computer is used. Monte Carlo method is a typical representative of the nonlinear global optimization method, but its weaknesses is low efficiency, which may require the order of thousands or even millions to obtain acceptable precision (Cheng et al., 2022b). With the rapid development of the current computer performance, the implementation of parallel computations enables some classic meta-heuristic algorithms, based on Monte Carlo ideas, such as simulated annealing, which is based on statistical mechanics, and the ant colony algorithm, which is based on biology (Stoffa and Sen, 1991; Mallick, 1995), in order for the inversion of actual seismic data. Ant colony optimization (ACO) is an iterative algorithm based on the natural system; it is inherently highly parallel because it simulates a range of solutions and has no strict restrictions on the initial model, but the efficiency of ACO is influenced by the population size heavily, that is, convergence can be premature if a small number of models are employed, and ACO does not always guarantee it will converge to the global optimal solution within a limited time.

Quantum computing methods have been hailed as the future of computing sciences (Moradi et al., 2018). This can be attributed to their theoretically proven supercomputing speed, better stability, and effectiveness. Although many studies have been conducted to examine the applicability of quantum computing in general (e.g., Moradi et al., 2018; Liu et al., 2018), its potential in geophysics has received limited attention. Quantum computing principles such as superposition, entanglement, and quantum tunneling can bring a paradigm shift in computing by achieving substantial speed-up over its classical counterpart (Hoos and Stϋtzle, 2018; Kallel et al., 2013). Quantum computing offers reversible computational logic, which helps in the minimization of power consumption while preserving computational information, making it a good alternative for reducing overall computational complexity and resource overhead. Additionally, quantum superposition using the Hadamard gate will allow parallel exploitation of the search space and thus may offer the desired solution faster than classical ACO synthesis through quantum interference (Grover, 1996). This improvement in QACO results in improved convergence efficiency and accuracy of results in optimization problems (Lahoz-Beltra, 2016). Despite these improvements, QACO shows two major deficiencies when solving relatively more complex problems: first, the slow rate of convergence. The angle of rotation, which is usually fixed, and the direction of the quantum rotating gate are obtained from a look-up table (Han et al., 2001), making the algorithm less flexible and resulting in a slower convergence rate in the early stages. Second, there is high tendency of being trapped in a local minimum. QACO uses the quantum rotating gate to update genes with the best fit. This increases the population of the best-fit genes but ultimately results in the loss of diversity, which increases the tendency of being trapped in a local minimum.

To address these shortcomings in the QACO method, in this study, we propose an improved version called the hybrid quantum ant colony optimization (HQACO). It combines a self-adaptive search strategy and the operations of the quantum gate, enjoying the advantages of quantum computing and ACO. This positions the HQACO as a great tool for global optimization using a small population size. It also enjoys the advantages of the results being independent of the initial model and high likelihood of obtaining the global solution. We verified its reliability and stability by conducting synthetic tests using models based on synthetic and actual logging data. The results show that HQACO is a strong candidate for global optimization with fast convergence speed and robust stability.

2 Methodology

2.1 Forward modeling

Based on the three-dimensional wave equation and the boundary conditions of continuous displacement and stress at the interface of the medium, Zoeppritz (1919) deduced the expression of the reflection and transmission coefficient when the wave propagates to the interface, which is described as the Zoeppritz equation. Assuming a solid–solid interface between two homogeneous isotropic elastic half-spaces, P-wave velocity, S-wave velocity, and density of the upper half-space are denoted by $V_{P1},V_{S1},{\rho }_1$, respectively; P-wave velocity, S-wave velocity, and density of the lower half-space are denoted by $V_{P2},V_{S2},{\rho }_2$, respectively; $i_1$ and $i_2$ are the incident angle and transmitted angle of P-wave; $j_1$ and $j_2$ are the incident angle and transmitted angle, respectively, of S-wave; the ray parameter is constant and described as $p=\sin i_1/V_{P1}=\sin i_2/V_{P2}=\sin j_1/V_{S1}=\sin j_2/V_{S2}$. Aki and Richards (2002) derived the accurate solutions of the Zoeppritz equation as

$R_{PP}=\frac{\left[\left(b\frac{\mathit{cos}{i}_1}{V_{P1}}-c\frac{\mathit{cos}{i}_2}{V_{P2}}\right)F-\left(a+d\frac{\mathit{cos}{i}_1}{V_{P1}}\frac{\mathit{cos}{j}_2}{V_{S2}}\right)H{p}^2\right]}{D},$

where

$\left\{\begin{array}{l}\begin{array}{l}a={\rho }_2\left(1-2V_{S2}^2{p}^2\right)-{\rho }_1\left(1-2V_{S1}^2{p}^2\right)\\ b={\rho }_2\left(1-2V_{S2}^2{p}^2\right)+2{\rho }_1{V}_{S1}^2{p}^2\end{array}\\ \begin{array}{l}c={\rho }_1\left(1-2V_{S1}^2{p}^2\right)+2{\rho }_2{V}_{S2}^2{p}^2\\ d=2\left({\rho }_2{V}_{S2}^2-{\rho }_1{V}_{S1}^2\right)\end{array}\end{array}\right.$

and

$\left\{\begin{array}{l}\begin{array}{c}E=b\frac{\mathit{cos}{i}_1}{V_{P1}}+c\frac{\mathit{cos}{i}_2}{V_{P2}}\\ F=b\frac{\mathit{cos}{j}_1}{V_{S1}}+c\frac{\mathit{cos}{j}_2}{V_{S2}}\\ G=a-d\frac{\mathit{cos}{i}_1}{V_{P1}}\frac{\mathit{cos}{j}_2}{V_{S2}}\end{array}\\ H=a-d\frac{\mathit{cos}{i}_2}{V_{P2}}\frac{\mathit{cos}{j}_1}{V_{S1}}\\ D=EF+GH{p}^2\end{array}\right..$

As shown above, although the Zoeppritz equation is a highly nonlinear function with respect to these properties, the formulations are explicitly valid in the isotropic media, and it uses few numbers of approximations, which allows the prediction of the reflection coefficients to be accurate from near to far incident angles. Then, global optimization schemes are suggested to be involved to treat AVO inversion. Combining the exact Zoeppritz equation with the global optimization methods can help inversion converge to global solutions.

2.2 Conventional objective optimization

Nonlinear inversions are usually done by minimizing an objective function. The L2 norm and cross-correlation coefficient are commonly objective functions to define an error, misfit, or similarity in the optimization problem. Then, we can obtain the objective function of the form as follows:

$m_{INV}=\mathit{arg}{\displaystyle \underset{m}{\min }}\left[{\left(r-G\left(m\right)\right)}^T\left(r-G\left(m\right)\right)\right],$

where r is the data vector, denoted as $r={\left[r_1,r_2,\dots ,r_n\right]}^T$, and model parameters $m={\left[m_1,m_2,\dots ,m_p\right]}^T$ are all positive; and $n$ and $p$ are the numbers of observed data points and model parameters, respectively. Arg min returns the value of m, which minimizes the given function. It should be noted that these objective functions are the simplest form and are affected by the energy of noise.

2.3 Hybrid quantum ant colony optimization

Simple ACO functions in two operating modes: forward (from the nest toward the food) and backward (from food back to the nest). Forward ants build a solution by probabilistically choosing the next node to move to among those in the adjacent positions with respect to the current node. This probabilistic choice is biased by pheromone trails previously deposited on the paths by other ants. Forward ants do not deposit any pheromone, which, when associated with deterministic backward moves, helps eliminate loop formation.

First, an ant k locating on node i moves to j as the next node using the following constraint for the transition function:

$m\left(j\right)=\left\{\begin{array}{l}\mathit{arg}\left\{\mathit{max}\left[{\left({\tau }_{ij}\left(t\right)\right)}^{\alpha }\bullet {\left({\eta }_{ij}\right)}^{\beta }\right]\right\},r\le r_0\\ j^{\prime },r>r_0\end{array}\right.,$

where ${\tau }_{ij}$ is the pheromone trail at time t, ${\eta }_{ij}$ is the problem-specific heuristic information, $\alpha$ and $\beta$ are impacts of pheromone trail and heuristic information, respectively, r is the random number with uniform distribution in [0, 1], r0 is the pre-specified parameter ranging from 0 to 1, and $j^{\prime }$ is the target point selected according to the following probability distribution:

$P_{ij}^k\left(t\right)=\frac{{\left({\tau }_{ij}\left(t\right)\right)}^{\alpha }\bullet {\left({\eta }_{ij}\right)}^{\beta }}{{\displaystyle \sum _{j\in N_i^k}}{\left({\tau }_{ij}\left(t\right)\right)}^{\alpha }\bullet {\left({\eta }_{ij}\right)}^{\beta }},$

where $N_i^k$ is the neighborhood of ant k when located on node i. Ant k traversing in the backward mode through the arc (i; j) will update the pheromone value as follows:

${\tau }_{ij}^{*}=\left(1-\rho \right){\tau }_{ij}+{\displaystyle \sum _{k=1}^K}{∆\tau }_{ij}^k,$

where ${∆\tau }_{ij}^k$ indicates the amount of pheromone deposition to all arcs, $\left.\rho \in (0,1\right]$ is the pheromone evaporation parameter, and K indicates the amount of ants.

In a quantum computer, a two-level quantum system called qubit acts as the physical media to store the information units. A qubit is the basic unit of information in quantum computation and is described by a superposition of the basis states $\phi ⟩⟩=\alpha \left|\left|0⟩\right.⟩\right.+\beta \left|\left|1⟩\right.⟩\right.$, where $\left|0⟩〉\right.$ and $\left|1⟩〉\right.$ are denoted as two basic states (Han and Kim, 2000); α and β are complex numbers and satisfy ${\left|\alpha \right|}^2+{\left|\beta \right|}^2=1$. ${\left|\alpha \right|}^2$ and ${\left|\beta \right|}^2$ are called the probability amplitude of the corresponding states $\left|0⟩〉\right.$ and $\left|1⟩〉\right.$ of qubit, respectively. We can note that a qubit is not a value of 0 or 1, but there is a possibility that it can be represented with values. A qubit can contain the information of both state |0⟩ and |1⟩, and a superposition state $\phi ⟩⟩$ also can be represented by a unit vector of a two-dimensional Hilbert space as $\left|\phi ⟩\right.〉=\left[\begin{array}{c}\alpha \\ \beta \end{array}\right]$ (Williams, 2010).

The quantum rotating gate is the core operator in the evolution operation, and it directly affects the performance of the algorithm. As model parameters are in superposition states, the genetic qubits in the population should be updated by the quantum rotating gates to adjust the probability amplitude and constitute new individuals. Quantum rotating gates can be designed according to the practical problems and can be usually defined as $G=\left[\begin{array}{cc}\mathit{cos}\theta & -\mathit{sin}\theta \\ \mathit{sin}\theta & \mathit{cos}\theta \end{array}\right]$, where θ is the rotating angle. The rotation strategy adopted is given by the following equation: $\left|{\phi }^{\prime }⟩\right.〉=G\left(\theta \right)\times \left|\phi ⟩\right.〉=\left[\begin{array}{cc}\mathit{cos}\theta & -\mathit{sin}\theta \\ \mathit{sin}\theta & \mathit{cos}\theta \end{array}\right]\left[\begin{array}{c}\alpha \\ \beta \end{array}\right]=\left[\begin{array}{c}{\alpha }^{\prime }\\ {\beta }^{\prime }\end{array}\right]$, where $\left|{\phi }^{\prime }⟩\right.〉$ is the updated quantum superposition state and α ‘and β' are the probability amplitudes of the quantum state after rotation. Generally, the rotation angle and rotation direction of the quantum gates are empirically determined in advance (Layeb and Saidouni, 2007; Han and Kim, 2000), and the adjustment strategy of the rotation angle and rotation direction for the quantum gate is updated in accordance with the self-adaptive rotating strategy proposed by Cheng et al. (2022b). A dynamic rotating method to adjust the rotation angle of the quantum rotating gate is proposed as follows:

$∆\theta =-sgn\left(A\right)\times \left({\theta }_{\mathit{max}}-{\theta }_{\mathit{min}}\right)\times \left|\frac{F_{Cur}-F_{\mathit{min}}}{C}\right|,$

where $A=\left|\begin{array}{cc}{\alpha }_B& {\alpha }_i\\ {\beta }_B& {\beta }_i\end{array}\right|$, ${\alpha }_B,{\beta }_B$ are the probability amplitudes of a qubit in the currently optimal solution, and ${\alpha }_i,{\beta }_i$ are the probability amplitudes of the corresponding qubits in the current search solution. Then, the direction of the rotation angle is taken as follows: when $A\ne 0$, the direction is $-sgn\left(A\right)$; when $A=0$, the direction can be positive or negative. ${\theta }_{\mathit{max}}$ represents the maximum value in the interval, ${\theta }_{\mathit{min}}$ represents the minimum value in the interval, $F_{\mathit{min}}$ represents the fitness of the best optimal individual, $F_{Cur}$ represents the fitness of the current individual, and C is a positive integer constant. The proposed method can associate the rotating angle with the fitness value and dynamically adjust the rotating angle, thereby improving the total convergence rate without losing the accuracy. Adjusting the corresponding qubit based on the fitness can make the individual update toward the direction of the optimal model parameter, and the rotating angle gradually reduces as the fitness converges to the optimal solution. Results listed in Table 1 represent the average evaluation number requested to the global minimum, which reduces from 321 to 217 after applying the self-adaptive rotating strategy, and the optimal solution is closer to the global minimum.

Table 1

Table 1. Average numbers of model evaluations required to converge to the global minimum and the optimal result after improvements.

Quantum gate, in practical applications, is an assistant operation with the purpose to avoid losing important information of the population and enhance the local searching ability, as the diversity maintained by the superposition state is not enough to reflect the superiority of the quantum algorithm and premature convergence may occur due to the effect of noise. Therefore, quantum gates will be added to improve the performance of the algorithm.

The quantum NOT gate (Pauli-X gate) is commonly adopted to implement quantum gates. Through exchanging the probability amplitude of each gene $\left|0⟩〉\right.$ and $\left|1⟩〉\right.$, the chance of the obtained result after measurement is exchanged. Moreover, the quantum NOT gate, as the quantum gate operator, rotates the angle of the gene by $\frac{\pi }{2}-2\theta$. As the rotation angle is slightly large, the convergence speed of the algorithm can be accelerated in the early stage of the quantum genetic algorithm, and at the same time, the diversity of the population can be guaranteed. However, the quantum NOT gate may cause premature phenomenon in the late stage. Consequently, we consider using Hadamard gates for quantum gates when the convergence is in the late stage. The matrix representation of the Hadamard gate is $H=\frac{1}{\sqrt{2}}\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]$, and the mutation of genes can be represented as

$H\left(\left|\mathrm{\psi }⟩|〉\right.\right)=\frac{1}{\sqrt{2}}\left[\begin{array}{cc}1& 1\\ 1& -1\end{array}\right]\left[\begin{array}{c}\mathit{cos}\theta \\ \mathit{sin}\theta \end{array}\right]=\left[\begin{array}{c}\mathit{cos}\left(\frac{\pi }{4}-\theta \right)\\ \mathit{sin}\left(\frac{\pi }{4}-\theta \right)\end{array}\right]=\left[\begin{array}{c}\mathit{cos}\left(\theta +\frac{\pi }{4}-2\theta \right)\\ \mathit{sin}\left(\theta +\frac{\pi }{4}-2\theta \right)\end{array}\right].$

From the formula, we can find that the angle of the qubit is rotated by $\frac{\pi }{4}-2\theta$. Therefore, the Hadamard gate is more suitable as a mutation operator when the individual is close to the optimal value, to prevent the information of the current optimal individual from disappearing and to ensure the stability of the population. Table 1 shows that using quantum gates leads to higher efficiency, the average calculation number decreases to 151, and the higher accuracy is when optimal solution closes to 0.

3 Inversion test

We present two numerical examples to examine the effectiveness of our proposed algorithm. The misfit between the predicted and observed data is used as the objective function to be minimized. An objective function with multiple local minima reflects high ill-conditioning of the inversion problem. Several performance criteria are used to evaluate the performance of the algorithms, such as the number of iterations and robustness. In all synthetic examples, the incidence angles for the AVA gathers are 1°–40°, and the wavelets used for the forward modeling are Ricker wavelets. First, we design a 1-dimensional 13-layer horizontal theoretical stratigraphic model to test the accuracy and efficiency of the 3 optimization algorithms: ACO, QACO, and HQACO under different search windows. Next, the actual logging data are used to generate synthetic seismic records, and the convergence efficiency and robustness of the inversion algorithm are tested by adding random noise.

3.1 Numerical test

3.1.1 Example I

First, we build a 13-layer horizontal layered model to test the stability and accuracy of HQACO, QACO, and ACO in different search windows. Given the termination criterion, maximum generation 500 and error tolerance 10–5, two kinds of search windows with three different search ranges, are used to detail the influence and choice of search windows required for inversions without any noise. Figure 1A shows the information of model I. Figure 1B shows the amplitude versus angle gathers generated by exact Zoeppritz forward modeling from 1° to 40°. Figure 2 shows the inversion result under the fixed search window with the deviation of 15%, 10%, and 5%, respectively. Figure 3 shows the result of the optimized search window with the deviation of 15%, 10%, and 5%. The optimized search window represents the upper and lower constraint bounds, which can be obtained from the initial models in practice. It is also noted that the dynamic search window in the nonlinear inverse problem is not applied to offer the low frequency component because the low frequency trend can be obtained in the static window. Instead, the dynamic window is used for a high computational efficiency and resolution. With the same length of chromosome, the narrow-bound constraint results in high resolution. Comparing the results shown in Figures 1, 2, it can be observed that the inversion result by HQACO has the highest accuracy among the three algorithms, and the accuracy of QACO is higher than that of ACO. In addition, the results obtained by the dynamic search window show better performance than those obtained by the static search window, and the inversion results increasingly deviate from the true models as the size of the search window expands. The relative errors (REs) of the results are calculated and shown in Tables 2–4, and all of them can conclude the consistent results with the figures. It is clear that HQACO is suitable for the inversion in such a theoretical situation, and the dynamic search window with narrow bounds can significantly decrease the inversion uncertainty. Then, we test the efficiency of convergence for different algorithms. We run each algorithm with the parameter range of 5%, and the maximum number of model evaluations is set to 1,000. Figure 4 shows that the results of all of the methods are close enough to the global minimum with limited iteration numbers, whereas ACO shows the lowest convergence rate because it has the highest computational complexity, which means minimizing computational cost may cause undesirable solutions.

Figure 1

Figure 1. (A) P-wave velocity, S-wave velocity, and density of model I are shown from left to right. (B) Synthetic seismic angle gathers of model I.

Figure 2

Figure 2. Using static search windows to constrain the AVO inversion with different bounds: (A) 15% deviation, (B) 10% deviation, and (C) 5% deviation. P-wave velocity, S-wave velocity, and density are shown from left to right. The solid black lines denote the true models. The red dashed curves denote the inverted models by HQACO, where the upper and lower bounds of search windows are denoted by gray dashed curves. The dashed blue lines and green lines indicate the inverted results by QACO and ACO, respectively.

Figure 3

Figure 3. Using dynamic search windows to constrain the AVO inversion with different bounds: (A) 15% deviation, (B) 10% deviation, and (C) 5% deviation. P-wave velocity, S-wave velocity, and density are shown from left to right. The solid black lines denote the true models. The red dashed curves denote the inverted models by HQACO, where the upper and lower bounds of search windows are denoted by gray dashed curves. The dashed blue lines and green lines indicate the inverted results by QACO and ACO, respectively.

Table 2

Table 2. Relative errors between the inversion results and true models by using ACO with different search windows.

Table 3

Table 3. REs between the inversion results and true models using QACO with different search windows.

Table 4

Table 4. REs between the inversion results and true models using HQACO with different search windows.

Figure 4

Figure 4. Evolution process with ACO (green), QACO (blue), and HQACO (red) using the 5% deviation dynamic search window.

3.1.2 Example II

Based on the well-log information, we calculate the PP reflection coefficients using the exact Zoeppritz equation in the time domain and then generate the PP AVA gathers by convolving the reflectivity with the real wavelets in different angles. The objective function is the L2 norm between the predicted and real angle gathers. The maximum number of fitness evaluations that we allow for algorithms to minimize the error is 500*N, where N is the dimension of the problem. In addition, different noise energies are added to the synthetic data to test the robustness of the proposed algorithm. Algorithm efficiency, like the number of function evaluations and the number of iterations, are likely to be correlated with the CPU time. Figure 5 shows the P-wave velocity, S-wave velocity, and density of well log. Figure 6A represents the AVA gathers computed using the Zoeppritz equation, and random noise with a signal-to-noise ratio (SNR) of 4 and 2 is added to the synthetic gathers, as shown in Figures 6B, C. Then, the results inverted by ACO and QACO with different SNR are shown in the Figure 7. It can be seen that the inversion results by QACO show higher accuracy than those using ACO; however, the error of density badly increases with the enhancement of noise energy. The results inverted using QACO and HQACO with different SNR are shown in Figure 8, and inversion results using HQACO without noise are closely consistent with real-model parameters. Moreover, results using HQACO show higher accuracy than those using QACO in different SNR. The relative errors (REs) and correlation coefficients (CCs) of the results are calculated and shown in Tables 5–7, and the lower SNR makes relative errors become larger, whereas the lower SNR makes correlation coefficients become smaller. Figures 9–11 show the residual between the predicted seismic data and the observed data (Figure 6) with different signal-to-noise ratios. It is shown that the difference between the observed model and predicted model is closely related to the magnitude of noise, where the error of HQACO is the smallest and the results using ACO show larger errors, especially at large offsets. Moreover, the convergence graph is a useful tool to show the average error performance of the total runs, and the average iteration curves shown in Figure 12 express that HQACO can spend less time converging to global minimum and reach the highest accuracy among the three algorithms, followed by QACO, and ACO is the worst. Above all, the proposed algorithm has higher accuracy, higher efficiency, and stronger robustness than the other two algorithms, and Bayesian inference is recommended to ensure the stability of inverted density.

Figure 5

Figure 5. Model II: well data of real data are used to generate the synthetic data for inversion. P-wave velocity, S-wave velocity, and density are shown from left to right.

Figure 6

Figure 6. Synthetic angle gathers with different levels of noise: (A) without noise, (B) signal-to-noise ratio = 4, and (C) signal-to-noise ratio = 2.

Figure 7

Figure 7. Comparison among the real logs (black solid), inversion results by QACO (blue), and inversion results by ACO (green) for synthetic model II with different levels of noise: (A) without noise, (B) signal-to-noise ratio = 4, and (C) signal-to-noise ratio = 2. P-wave velocity, S-wave velocity, and density are shown from left to right. The upper and lower bounds of search windows are denoted by gray curves.

Figure 8

Figure 8. Comparison among the real logs (black solid), inversion results by HQACO (red), and inversion results by QACO (blue) for synthetic model II with different levels of noise: (A) without noise, (B) signal-to-noise ratio = 4, and (C) signal-to-noise ratio = 2. P-wave velocity, S-wave velocity, and density are shown from left to right. The upper and lower bounds of search windows are denoted by gray curves.

Table 5

Table 5. REs and CCs between the inversion results and true models using ACO with different SNR.

Table 6

Table 6. REs and CCs between the inversion results and true models using QACO with different SNR.

Table 7

Table 7. REs and CCs between the inversion results and true models using HQACO with different SNR.

Figure 9

Figure 9. Residual between input seismic data and predicted seismograms with different methods when without noise: (A) ACO, (B) QACO, and (C) HQACO.

Figure 10

Figure 10. Residual between input seismic data and predicted seismograms with different methods when signal-to-noise ratio = 4: (A) ACO, (B) QACO, and (C) HQACO.

Figure 11

Figure 11. Residual between input seismic data and predicted seismograms with different methods when signal-to-noise ratio = 2: (A) ACO, (B) QACO, and (C) HQACO.

Figure 12

Figure 12. Evolution process with ACO (green), QACO (blue), and HQACO (red).

4 Discussion and conclusion

We propose a nonlinear AVO inversion method based on the exact Zoeppritz equation and apply it to the synthetic model. When dealing with high-contrast interfaces with far-offset information, the new inversion method based on the accurate Zoeppritz equation is sufficient to achieve the promising results. QACO shows greater global search competence and higher efficiency than ACO because the quantum qubit increases the population diversity and the quantum rotating gate is involved. HQACO, combining the quantum rotating gate with the self-adaptive rotating strategy and quantum gate, can adaptively update the search step and accelerate the convergence. All inversions use a dynamic search window, which can be constructed by combining the interpreted horizon with smoothed logs calibrated using seismic data. Inversion using dynamic search windows is more efficient than that using static search windows. Although it is recommended that special care must be taken in setting the search window range, in general, larger search windows require a larger computation time and have a greater likelihood of obtaining a global solution, whereas smaller search windows may miss the true solution. In conclusion, although the nonlinear inversion method based on the exact equations is time-consuming, HQACO is reliable in accuracy for reservoirs with strong contrast, long offset ranges, and noise energy using the exact Zoeppritz equation. The proposed method also significantly improves the inversion efficiency. The current quantum algorithm is developed in the current computer, whereas quantum computers promise computational improvement for a wide variety of applications, indicating that the algorithm should be designed for quantum computers to make possible further development of the efficiency of the inversion method. The proposed inversion method in conceptual models has been tested, and it can be applied in various geophysics problems and disciplines. With the recent breakthroughs in the construction of the universal quantum computer, future work will be devoted to incorporating multicomponent information, anisotropy, and complex wave propagation effects.

Data availability statement

The original contributions presented in the study are included in the article; further inquiries can be directed to the corresponding author.

Author contributions

YC: conceptualization, data curation, formal analysis, funding acquisition, investigation, methodology, project administration, resources, software, supervision, validation, visualization, writing–original draft, and writing–review and editing. JC: conceptualization, supervision, and writing–review and editing. HM: supervision and writing–review and editing. DF: writing–review and editing. CW: validation and writing–review and editing. CZ: writing–review and editing. MS: supervision and writing–review and editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. The study was sponsored by the National Natural Science Foundation of China (42304134).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Aki, K., and Richards, P. G., 2002, Quantitative seismology.

Google Scholar

Aleardi, M., and Mazzotti, A. (2017). 1D elastic full-waveform inversion and uncertainty estimation by means of a hybrid genetic algorithm–Gibbs sampler approach. Geophys. Prospect. 65 (1), 64–85. doi:10.1111/1365-2478.12397

CrossRef Full Text | Google Scholar

Aleardi, M., Pierini, S., and Sajeva, A. (2019). Assessing the performances of recent global search algorithms using analytic objective functions and seismic optimization problems. Geophysics 84 (5), R767–R781. doi:10.1190/geo2019-0111.1

CrossRef Full Text | Google Scholar

Buland, A., and Omre, H. (2003). Bayesian linearized AVO inversion. Geophysics 68 (1), 185–198. doi:10.1190/1.1543206

CrossRef Full Text | Google Scholar

Cheng, J. W., Zhang, F., and Li, X. Y. (2022a). Seismic amplitude inversion for orthorhombic media based on a modified reflection coefficient approximation. Surv. Geophys. 43 (5), 1395–1433. doi:10.1007/s10712-022-09718-z

CrossRef Full Text | Google Scholar

Cheng, J. W., Zhang, F., and Li, X. Y. (2022b). Nonlinear amplitude inversion using a hybrid quantum genetic algorithm and the exact zoeppritz equation. Petroleum Sci. 19 (3), 1048–1064. doi:10.1016/j.petsci.2021.12.014

CrossRef Full Text | Google Scholar

Grover, L. K. (1996). “A fast quantum mechanical algorithm for database search,” in In Proceedings of the twenty-eighth annual ACM symposium on Theory of computing, Pennsylvania, Philadelphia, USA, May 22–24, 1996.

CrossRef Full Text | Google Scholar

Han, K. H., and Kim, J. H. (2000). “Genetic quantum algorithm and its application to combinatorial optimization problem,” in Proceedings of the 2000 congress on evolutionary computation. CEC00 (Cat. No. 00TH8512), La Jolla, CA, USA, 16–19 July, 2000 (IEEE), 2, 1354–1360.

Google Scholar

Han, K. H., Park, K. H., Lee, C. H., and Kim, J. H. (2001). “Parallel quantum-inspired genetic algorithm for combinatorial optimization problem,” in Proceedings of the 2001 congress on evolutionary computation (IEEE cat. No. 01TH8546), La Jolla, CA, USA, 16–19 July, 2000 (IEEE), 2, 1422–1429.

Google Scholar

Hoos, H. H., and Stϋtzle, T. (2018). “Stochastic local search,” in Handbook of approximation algorithms and metaheuristics (Chapman and Hall/CRC), 297–307.

CrossRef Full Text | Google Scholar

Kallel, L., Naudts, B., and Rogers, A. (2013). Theoretical aspects of evolutionary computing (Springer Science & Business Media).

Google Scholar

Lahoz-Beltra, R. (2016). Quantum genetic algorithms for computer scientists. Computers 5 (4), 24. doi:10.3390/computers5040024

CrossRef Full Text | Google Scholar

Layeb, A., and Saidouni, D. E. (2007). Quantum genetic algorithm for binary decision diagram ordering problem. Int. J. Comput. Sci. Netw. Secur. 7 (9), 130–135.

Google Scholar

Li, T., and Mallick, S. (2015). Multicomponent, multi-azimuth pre-stack seismic waveform inversion for azimuthally anisotropic media using a parallel and computationally efficient non-dominated sorting genetic algorithm. Geophys. J. Int. 200 (2), 1136–1154. doi:10.1093/gji/ggu445

CrossRef Full Text | Google Scholar

Liu, X., Cheng, J., Cai, Y., Mo, Q., Li, C., and Zu, S. (2022). Stochastic simulation of facies using deep convolutional generative adversarial network and image quilting. Mar. Petroleum Geol. 146, 105932. doi:10.1016/j.marpetgeo.2022.105932

CrossRef Full Text | Google Scholar

Liu, X., Chen, X., Cheng, J., Zhou, L., Chen, L., Li, C., et al. (2023). Simulation of complex geological architectures based on multi-stage generative adversarial networks integrating with attention mechanism and spectral normalizatio. IEEE Trans. Geosci. Remote Sens. doi:10.1109/TGRS.2023.3294493

CrossRef Full Text | Google Scholar

Liu, X., Li, J., Chen, X., Guo, K., Li, C., Zhou, L., et al. (2018). Stochastic inversion of facies and reservoir properties based on multi-point geostatistics. J. Geophys. Eng. 15 (6), 2455–2468. doi:10.1088/1742-2140/aac694

CrossRef Full Text | Google Scholar

Lu, J., Yang, Z., Wang, Y., and Shi, Y. (2015). Joint PP and PS AVA seismic inversion using exact Zoeppritz equations. Geophysics 80 (5), R239–R250. doi:10.1190/geo2014-0490.1

CrossRef Full Text | Google Scholar

Mallick, S. (1995). Model-based inversion of amplitude-variations-with-offset data using a genetic algorithm. Geophysics 60 (4), 939–954. doi:10.1190/1.1443860

CrossRef Full Text | Google Scholar

Mallick, S., and Adhikari, S. (2015). Amplitude-variation-with-offset and prestack-waveform inversion: a direct comparison using a real data example from the Rock Springs Uplift, Wyoming, USA. Geophysics 80 (2), B45–B59. doi:10.1190/geo2014-0233.1

CrossRef Full Text | Google Scholar

Moradi, S., Trad, D., and Innanen, K. A. (2018). “Quantum computing in geophysics: algorithms, computational costs, and future applications,” in SEG technical program expanded abstracts 2018 (Society of Exploration Geophysicists), 4649–4653. doi:10.1190/segam2018-2998507.1

CrossRef Full Text | Google Scholar

Ostrander, W. J. T. (1984). Plane-wave reflection coefficients for gas sands at nonnormal angles of incidence. Geophysics 49 (10), 1637–1648. doi:10.1190/1.1441571

CrossRef Full Text | Google Scholar

Russell, B. H., Gray, D., and Hampson, D. P. (2011). Linearized AVO and poroelasticity. Geophys. 76 (3), C19–C29. doi:10.1190/1.3555082

CrossRef Full Text | Google Scholar

Sen, M. K., and Stoffa, P. L. (2013). Global optimization methods in geophysical inversion. Cambridge University Press.

Google Scholar

Shuey, R. T. (1985). A simplification of the Zoeppritz equations. Geophysics 50 (4), 609–614. doi:10.1190/1.1441936

CrossRef Full Text | Google Scholar

Smith, G. C., and Gidlow, P. M. (1987). Weighted stacking for rock property estimation and detection of gas. Geophys. Prospect. 35 (9), 993–1014. doi:10.1111/j.1365-2478.1987.tb00856.x

CrossRef Full Text | Google Scholar

Stoffa, P. L., and Sen, M. K. (1991). Nonlinear multiparameter optimization using genetic algorithms: inversion of plane-wave seismograms. Geophysics 56 (11), 1794–1810. doi:10.1190/1.1442992

CrossRef Full Text | Google Scholar

Williams, C. P. (2010). Explorations in quantum computing. Springer Science & Business Media.

Google Scholar

Zhang, F., and Li, X. (2013). Generalized approximations of reflection coefficients in orthorhombic media. J. Geophys. Eng. 10 (5), 054004. doi:10.1088/1742-2132/10/5/054004

CrossRef Full Text | Google Scholar

Zhang, F., Zhang, T., and Li, X. Y. (2019). Seismic amplitude inversion for the transversely isotropic media with vertical axis of symmetry. Geophys. Prospect. 67 (9-Advances in Seismic Anisotropy), 2368–2385. doi:10.1111/1365-2478.12842

CrossRef Full Text | Google Scholar

Zhou, L., Li, J., Chen, X., Liu, X., and Chen, L. (2017). Prestack AVA inversion of exact Zoeppritz equations based on modified Trivariate Cauchy distribution. J. Appl. Geophys. 138, 80–90. doi:10.1016/j.jappgeo.2017.01.009

CrossRef Full Text | Google Scholar

Zoeppritz, K. (1919). On reflection and propagation of seismic waves. Gottinger Nachrichten, 1: 66–84.

Google Scholar

Zong, Z., Yin, X., and Wu, G. (2013). Elastic impedance parameterization and inversion with Young’s modulus and Poisson’s ratio. Geophysics 78 (6), N35–N42. doi:10.1190/geo2012-0529.1

CrossRef Full Text | Google Scholar