CN102854528B - Pre-stack non-linear inversion method based on particle swarm optimization algorithm - Google Patents

Abstract

The invention relates to a pre-stack non-linear inversion method based on a particle swarm optimization algorithm. Along with the sustained and steady growth of the economy of the country, the need for energy sources is also sustainedly and quickly increased and a serious test of shortage of oil and natural gas resources can inevitably appear. Reservoir prediction is an important research link in oil exploration production; information (contained in a seismic data) about structure, lithology and the like of subsurface formation needs to be fully utilized to find high-quality reservoirs; the information is frequently obtained through seismic inversion; but the conventional inversion is based on a linear and simplified model and is disturbed by a plurality of human factors, and high-precision reservoir parameters frequently cannot be obtained. According to the method, a research on the topic of the particle swarm optimization algorithm is carried out; the Zoeppritz equation is deduced and simplified on the basis of the idea and the principle of the particle swarm optimization algorithm to ensure that the equation is suitable for the solution of the particle swarm optimization algorithm; the particle swarm optimization algorithm is used for the pre-stack inversion; and the pre-stack inversion is respectively carried out on the two-dimensional theoretical model and a offshore seismic data so as to obtain a satisfactory result.

Description

Particle Swarm Optimization Algorithm Prestack Nonlinear Inversion Method

technical field

The invention relates to the technical field of oil and gas field exploration, belongs to the category of seismic data inversion, and specifically solves the non-linear seismic inversion problem through particle swarm optimization algorithm. the

Background technique

Reservoir lithology identification and fluid prediction using seismic data has always been the goal pursued by geophysicists. The theoretical basis of seismic prestack inversion (AVO inversion) is the famous Zoeppritz equations. Ostrander (1984) found that "the reflection amplitude of gas-bearing sandstone increases with the increase of offset, and the reflection amplitude of water-bearing sandstone decreases with the increase of offset" in the process of studying the seismic amplitude characteristics of "bright spot" sandstone reservoirs. This phenomenon has greatly improved the ability of hydrocarbon detection, drawn people's attention from post-stack to pre-stack, and marked the emergence of practical AVO technology. However, the Zoplitz equation is very complicated, and its physical meaning is not clear. Geophysicists at home and abroad have studied and simplified it in many forms, and obtained a series of approximate formulas. Bortfeld (1961) calculated the reflection coefficients of plane longitudinal waves and transmitted waves by using the method of approaching a single interface with the formation thickness approaching zero, and gave a simplified formula for distinguishing fluids and solids. The approximation equation proposed by Aki and Richards in 1980 focused on describing the influence of changes in the velocity and density of compressional and shear waves on the reflection coefficient. Shuey's (1985) approximation equation expresses the reflection coefficient as the relationship between the normal incidence and the sum of the near, middle and far different incidences, which more intuitively reflects the relationship between the amplitude and the incidence angle. Smith and Gidlow (1987) proposed a weighted superposition processing method for AVO analysis and introduced the concept of "fluid factor" when performing AVO analysis on gas-bearing sandstone. Starting from the Aki and Richards equations, Fatti et al. (1994) proposed a relatively simplified formula to predict the AVO response of gas-bearing sandstones under the premise of weakening the formation density term. The elastic impedance (EI) concept and mathematical expression proposed by Connolly (1999) greatly enriched the connotation of AVO analysis technology, expanded the extension of AVO inversion technology, and extended the traditional P-wave impedance inversion to sub-angle stacked data. Ozdemir et al. (2001) proposed the synchronous inversion of elastic parameters using AVO information, which is called "joint inversion". Ma Jinfeng et al. (2003) proposed the concept of generalized elastic wave impedance. In recent years, AVO analysis and stochastic inversion techniques have been continuously developed, pre-stack synchronous inversion methods have become increasingly mature, and non-linear inversion techniques are mostly used in synchronous inversion methods. AVO was originally proposed only to improve the ability of hydrocarbon detection. Today, the development of AVO has already gone beyond this category and penetrated into various fields of seismic exploration. It has been widely used in fracture detection, pressure prediction, reservoir dynamic detection, oil and gas prediction, and description of reservoir heterogeneity. the

Particle Swarm Optimization (PSO), also known as PSO, is a kind of Swarm Intelligence (SI). It was first proposed by American scholars Kennedy and Eberhart in 1995. Particle swarm optimization algorithm is inspired by the search strategy in the defense and predation behavior of bird groups. Many creatures in nature have certain group behaviors, such as flocks of birds, schools of fish, etc. Although a single individual in the group has only simple behavior rules, the behavior of the formed group is very complicated. Many scientists have conducted research, including computer simulations, on the group behavior of flocks of birds or fish. The particle swarm optimization algorithm is inspired from this model and used to solve the optimization problem, and the solution of the optimization problem is regarded as a point in the search space, which is called "particle". The movement of each particle is optimized according to the "flight experience" of itself and other particles, so as to achieve the purpose of searching the optimal solution in the whole space. Many scholars and researchers have proposed many improved PSO algorithms on the basis of the basic PSO algorithm in terms of parameter selection, topology structure and integration with other optimization algorithms. Kennedy and Eberhart proposed the Binary Particle Swarm Optimization (BPSO) algorithm in 1997, which is a discrete particle swarm (Discrete PSO, DPSO) algorithm based on continuous space. Clerc (2000) proposed the TSP-DPSO algorithm for the Traveling Salesman Problem (TSP), which is a DPSO based on discrete space. Jun Sun et al. (2003) introduced quantum behavior into particle swarm algorithm and proposed quantum particle swarm algorithm (Quantum PSO, QPSO). In 2004, Gao Ying proposed a particle swarm algorithm based on simulated annealing (SA-PSO). In 2007, P.S. Shelokar et al. proposed a hybrid algorithm based on ant colony and particle swarm optimization, that is, PSACO (Particle Swarm Ant Colony Optimization) algorithm. Zhou Yalan et al. (2008) introduced the idea of distribution estimation algorithm into the PSO algorithm, and proposed a discrete particle swarm (Estimation of Distribution PSO, EDPSO) algorithm based on distribution estimation. Lin Dongyi et al. (2008) introduced the immune mechanism into the BPSO algorithm, based on a certain attribute importance measure of the difference matrix of the decision table as the vaccine model, and proposed a minimum attribute reduction algorithm based on immune particle swarm optimization, called immune discrete Particle Swarm Optimization (IPSO Algorithm). The particle swarm optimization algorithm has developed rapidly in recent years and has been successfully applied to many fields such as function optimization, neural network training, pattern recognition and classification, fuzzy system control and engineering, and a large number of practical applications have proved its effectiveness. It has a good development prospect and deserves further research. the

Contents of the invention

The invention utilizes amplitude-preserving high-resolution seismic data to carry out nonlinear inversion, uses updated reflection coefficient formulas in the process of solving nonlinear problems, and uses particle swarm algorithm to search for optimization. This measure can reduce human error and effectively improve the accuracy of inversion. the

First make the following definitions:

In geophysics, the reflection coefficient of a continuous function x(t) is defined as,

$\underset{\mathrm{<span\; class="notranslate">\; \&Delta;t</span>}<span\; class="notranslate">\; \&Right\; Arrow;</span>\mathrm{<span\; class="notranslate">\; 0</span>}}{{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; \&Delta;x</span>}}{\mathrm{<span\; class="notranslate">\; 2</span>}\stackrel{<span\; class="notranslate">\; \&OverBar;</span>}{\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}}$

The continuous function x(t) is sampled into a discrete function at the sampling interval of dt, then,

${\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; x</span>}}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}{{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}$

Where x can be the longitudinal wave velocity α, the shear wave velocity β, and the density ρ. Subscript 1 indicates the parameters of the underlying formation, and subscript 2 indicates the parameters of the overlying formation. Defining the S-to-P wave velocity ratio

$\mathrm{<span\; class="notranslate">\; K</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}$

This patent uses the common reflection point gather in the reflection angle domain obtained by amplitude-preserving processing as input, and the technical scheme adopted to achieve the above purpose is as follows: 

Step 1: For any location underground. Firstly, the value ranges of P-wave velocity reflection coefficient, S-wave velocity reflection coefficient, density reflection coefficient, and S-wave velocity ratio are given, and the maximum number of iterations and population size are set. Given the inertia weight and learning factor, set the allowable error range. the

The inertia weight, learning factor, and population size have the following meanings:

Inertial weight: The balance between exploration ability and development ability is an important aspect that affects the performance of the optimization algorithm. For the particle swarm optimization algorithm, the balance of these two capabilities is achieved by the inertia weight w. A larger inertia weight makes the particle have a greater speed in its original direction, so it can fly farther in the original direction and has better exploration capabilities; a smaller inertia weight makes the particle inherit less speed in the original direction , so that the flight is closer and has better development capabilities. the

Learning factor: The learning factor is a non-negative constant, which represents the weight of particle preference, so that particles have the ability to self-summarize and learn from outstanding individuals in the group, so as to approach the optimal point in the group or in the neighborhood. Kennedy believes that the sum of the two learning factors should be about 4.0, and the search effect at this time is better. A common practice is to set them all to 2.05. the

Group size: When the group size is set to be small, the possibility of falling into a local optimum is high; when the group size is 1, the particle swarm optimization algorithm becomes a method based on individual search, and once trapped in a local optimum, it is impossible to jump out; When the group size is large, the optimization ability of the particle swarm optimization algorithm is better, but it will lead to a substantial increase in computing time, and when the number of groups grows to a certain level, continuing to grow will no longer have a significant effect. As for how many particles can participate in the search to achieve the desired effect, the predecessors calculated the average fitness for different populations by using multiple benchmark functions, and believed that the search efficiency is better when the population number is kept at about 30. the

Step 2: Initialize the P-wave velocity reflection coefficient, S-wave velocity reflection coefficient, Density reflection coefficient, and S-wave velocity ratio of each particle with random numbers, and set the initial velocity to zero, so that the individual optimal value and initial value of the particle The same, and calculate the fitness of each particle, choose the particle with the smallest fitness as the current global optimal solution. the

The calculation method of particle fitness is:

First define the matrix

${\mathrm{<span\; class="notranslate">\; m</span>}}^{<span\; class="notranslate">\; +</span>}<span\; class="notranslate">\; =</span>\left[\begin{array}{cccc}\mathrm{<span\; class="notranslate">\; D.</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; AD</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}\\ \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; KD</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; AD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="KBD"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">KBD</figure-callout></span>}\\ \mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \frac{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; K</span>}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="CB"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">CB</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; AD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \frac{\mathrm{<span\; class="notranslate">\; BC</span>}<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; ]</span>}{\mathrm{<span\; class="notranslate">\; K</span>}}\\ \mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="K"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">K</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; AC</span>}<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; ]</span>& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="CB"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">CB</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="DK"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">DK</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}\end{array}\right]$

A, B, C, D and K in the formula are defined in the following way

$\mathrm{<span\; class="notranslate">\; A</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}}{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}}$

$\mathrm{<span\; class="notranslate">\; B</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}}{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}}$

$\mathrm{<span\; class="notranslate">\; C</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}}{{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}}$

D = sini ₁

$\mathrm{<span\; class="notranslate">\; K</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&beta;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}{{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}$

where i ₁ is the angle of incidence. redefine

${\mathrm{<span\; class="notranslate">\; m</span>}}^{<span\; class="notranslate">\; -</span>}<span\; class="notranslate">\; =</span>\left[\begin{array}{cccc}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; D.</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; AD</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}\\ \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; D.</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; KD</span>}& \sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; AD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="KBD"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">KBD</figure-callout></span>}\\ \mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; D.</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \frac{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; K</span>}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="CB"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">CB</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; AD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& \frac{\mathrm{<span\; class="notranslate">\; BC</span>}<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; ]</span>}{\mathrm{<span\; class="notranslate">\; K</span>}}\\ \mathrm{<span\; class="notranslate">\; 2</span>}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="K"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">K</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\mathrm{<span\; class="notranslate">\; D.</span>}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; AC</span>}<span\; class="notranslate">\; [</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; ]</span>& <span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="CB"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">CB</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; <figure-callout\; id="2"\; label="DK"\; filenames="HSA00000748356000012.png"\; state="\{\{state\}\}">DK</figure-callout></span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}\sqrt{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; KBD</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}\end{array}\right]$

Then the solution of the Zoplitz equation can be obtained in the following form

${\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; pp</span>}}<span\; class="notranslate">\; =</span>\frac{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; m</span>}}^{<span\; class="notranslate">\; -</span>}<span\; class="notranslate">\; |</span>}{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; m</span>}}^{<span\; class="notranslate">\; +</span>}<span\; class="notranslate">\; |</span>}$

From this expression form, the expression of particle fitness can be obtained as

E＝||w(θ)*R _pp (R _α , R _β , R _ρ , K, θ)-d(θ)||

In the formula, E is particle fitness, w(θ) is seismic wavelet, and d(θ) is seismic observation data. the

Step 3: Compare whether the current global optimal fitness has reached the error tolerance range. If the error tolerance range is reached, the calculation is stopped, and the global optimum is output as the inversion result. Otherwise, go to step 4. the

Step 4: Update the particle's velocity, and update the particle's position with the new particle velocity. The number of iterations is increased once, and if the number of iterations exceeds the maximum number of iterations, the calculation is stopped and the global optimum is output. the

The update method of particle velocity and position is:

v _id (t+1)＝wv _id (t)+c ₁ r ₁ (p _id (t)-x _id (t))+c ₂ r ₂ (g _d (t)-x _id (t))

x _id (t+1)＝x _id (t)+v _id (t+1)

In the formula, v is the velocity and x is the position of the particle. m is the population size, n is the dimension of particles, 1≤i≤m, 1≤d≤n, w is called the inertia weight factor; c ₁ and c ₂ are called learning factors or acceleration factors, among which c ₁ is the adjustment The step size of the particle flying to its own best position, c ₂ is the step size for adjusting the particle to fly to the global best position; r ₁ and r ₂ are random numbers in [0, 1]; t is the current iteration algebra.

Step 5: Recalculate the fitness of the particles, and check whether the fitness of each particle is smaller than the fitness before updating. If it is less than, the individual optimum of the updated particle is the updated position of the particle. the

Step 6: Compare the new individual optimum, and select the individual optimum with the smallest fitness as the global optimum. Determine whether the global optimum reaches the tolerance range of the error, and output the global optimum as the inversion result if it reaches the tolerance range of the error, otherwise return to step 4. Until all calculation points underground are processed. the

The above specific embodiments are only used to illustrate the present invention, but not to limit the present invention. the

Description of drawings

Figure 1 is a schematic diagram of particle position update by particle swarm optimization algorithm. the

Figure 2 is a composite common reflection point gather in the reflection angle domain. the

Fig. 3 is a comparison of the P-wave reflection coefficient obtained by inversion of the gather shown in Fig. 2 using this method and the model P-wave reflection coefficient. the

Fig. 4(a) is the comparison of the P-wave velocity obtained by inversion of the gather shown in Fig. 2 with the model P-wave velocity;

Figure 4(b) is the comparison of the shear wave velocity obtained by inversion of the gather shown in Figure 2 with the model shear wave velocity;

Figure 4(c) is the comparison between the density obtained by inversion of the gather shown in Figure 2 and the model density;

Fig. 4(d) is a comparison between the S-P wave velocity ratio obtained by the inversion of the gather shown in Fig. 2 and the model S-P wave velocity ratio. the

Figure 5 is a comparison of the P-to-S wave velocity ratio (right) obtained by inversion of the method in the actual data in a certain area and the P-to-S wave velocity ratio (left) obtained by commercial software. the

Detailed ways

Illustrate with a synthetic gather:

Firstly, the Zoplitz equation was used to perform forward modeling on the model, and the reflection coefficients of each layer within the range of 1° to 30° were calculated every 3°. The calculated 10 reflection coefficient sequences are respectively convoluted with the 40Hz zero-phase Reker wavelet, and the obtained synthetic seismic records, that is, angle gather data, are shown in Fig. 2 . the

These ten angle gather data are used as input for inversion. The minimum error function E between the calculated P-wave reflection coefficient R _pp and the actual seismic record is taken as the objective function. The steps of pre-stack inversion using particle swarm optimization algorithm are as follows:

(1) The angle gather record obtained by the forward modeling of the model is taken as the observed value (true value). the

(2) Use the particle swarm optimization algorithm to generate the reflection coefficient sequence of each angle gather, such as r(t) candidate solutions. The position vector of each particle is the candidate solution of r(t), and the dimension of each particle is 4, that is, the four independent variables, longitudinal wave reflection coefficient R _α , shear wave reflection coefficient R _β , density reflection coefficient R _ρ , and transverse and longitudinal wave Speed than K.

(3) Calculate the fitness function E according to step 2. Comparing the calculated R _pp with the angle gather data, calculate the fitness function E.

(4) Substitute the obtained fitness function E into the global optimal solution G _best , and update the position of each particle according to step 4.

(5) Use the particle swarm optimization algorithm to iterate cyclically to continuously optimize the optimal solution until the termination condition is met (such as E is less than a certain truncation error, or reaches the maximum number of iterations), then the individual optimum is consistent with the global optimum at this time, and all The particles converge to one point, and the corresponding position vector is the required longitudinal wave reflection coefficient R _α , shear wave reflection coefficient R _β , density reflection coefficient R _ρ , and shear-to-pitch wave velocity ratio K.

In this example, when the particle swarm optimization algorithm is used to invert the reflection coefficient of the corner gather data of the numerical model, the parameter values are as follows: the error function is less than 1e-6, the population size m=30, the inertia factor w=0.1, and the learning factor c ₁ =c ₂ = 2.04, the maximum value of the particle position is [0.1, 0.1, 0.1, 0.75], and the minimum value is [-0.1, -0.1, -0.1, 0.35], corresponding to the longitudinal wave reflection coefficient R _α and the shear wave reflection coefficient R _β , the density reflection coefficient R _ρ , and the limit range of the shear-to-pillar velocity ratio K.

The comparison between the inverted reflection coefficient and the original reflection coefficient is shown in Figure 3. It can be seen that the inverted reflection coefficient (dashed line) agrees well with the true reflection coefficient (circle). Through the reflection coefficients R _α , R _β , R _ρ , and K obtained from the inversion, elastic parameters (solid lines) such as P-wave velocity, S-wave velocity, and density of each layer can be calculated. The comparison with the elastic parameters (dotted line) of the model is shown in Fig. 4. It can be seen from Figure 4 that the P-wave data and density data are completely consistent with the real situation, but there is a very small error between the S-wave data and the P-S wave velocity ratio. The reason for these errors is that in the inversion process, the accuracy of the S-wave component Affected by factors such as initial value setting and petrophysical relations, the stability is not as good as that of the longitudinal wave component; the second is due to the truncation error of the computer during the inversion process. In order to reduce the error as much as possible, not only the longitudinal wave data are used as the inversion observation value, but also the converted shear wave information can be used when possible. Even in the case of having VSP data, the data of reflected longitudinal wave, reflected shear wave, transmitted longitudinal wave, and transmitted shear wave can be used as inversion observations, then the various reflection coefficients obtained from the inversion and the rock’s elasticity derived from it parameters will be more accurate.

references

Zoeppritz K, Erdbenwellen VIIIB. On the reflection and propagation of seismic waves. Gottinger Naehriehten, 1919, 66-84

Ostrander W J. Plane waves reflection coefficients for gas sands at normal angles of incidence. Geophysics, 1984, 49(10): 1637-1648.

Bortfeld R. Approximation to the reflection and transmission coefficients of plane longitudinal and transverse wave. Geophysics Prospecting, 1961, 9(4): 485-502.

Aki K, Richards P G. Quantitative seismology: Theory and methods. W.H. Freeman and Co, 1980, 100-170.

Shuey. A simplification of the Zoeppritz equations. Geophysics, 1985, 50(4): 609-634.

Smith G C, Gidlow P M. Weighted stacking for rock property estimation and detection of gas. Geophysical Prospecting, 1987, 35: 993-1014.

Fatti J L, Vail P J, Smith G C, et al. Detection of gas in sandstone reservoirs using AVO analysis: A case seismic case history using the Geostack technique. Geophysics, 1994, 59(5): 1362～1376.

Connolly P. Elastic Impedance. The Leading Edge, 1999, 18(4): 438～452.

Ostrander W J. Plane-wave reflection coefficients for gas sands at non-normal angels of incidence. Exploration Geophysics, 1984, 49: 1637-1648.

Ma Jinfeng. Forward and inversion of generalized elastic impedance in seismic exploration. Acta Geophysics, 2003, 46(1): 118-124.

Du Yan. AVO technology and its application in gas-bearing detection of Sulige gas field reservoirs: (Master's thesis). Chengdu. Chengdu University of Technology, 2009.

Kennedy, Eberhart R C. Particle Swarm Optimization. International conference on Neural Networks, 1995, 1942-1948.

Liu Bo. Particle Swarm Optimization Algorithm and Its Engineering Application. Beijing: Electronic Industry Press, 2010. 

Kennedy J, Eberhart R C.A discrete binary version of the Particle Swarm Algorithm. The 1997 Conference on System, Cybernetics and Informatics, 1997, 4104-4109.

Clerc M. Discrete Particle Swarm Optimization, Illustrated by the Traveling Salesman Problem. http://www.mau-ricecierc.net, 2000. 

Sun J, Feng B, Xu W. Particle Swarm Optimization with Particle Having Quantum Behavior. IEEE Proc.Of Congress on Evolutionar Computation, 2004, 325-331.

Gao Ying, Xie Shengli. Particle Swarm Optimization Algorithm Based on Simulated Annealing. Computer Engineering and Application, 2004(1): 47-50.

Shelokar P S, Siarry P, Jayaraman V K, Kulkarni B D. Particle swarm and ant colony algorithms hybridized for improved continuous optimization. Applied Mathematics and Computation, 2007, 188: 129-142.

Zhou Yalan, Wang Jiahai, etc. A Discrete Particle Swarm Optimization Algorithm Based on Distribution Estimation. Acta Electronica Sinica, 2008, 36(6): 1242-1248.

Lin Dongyi, Liao Jiankun. An Immune Discrete Particle Swarm Optimization Algorithm for the Minimum Reduction Problem. Small Microcomputer Systems, 2008, 29(6): 1089-1092.

He Yi. Research on Nonlinear Inversion and Prediction of Seismic Reservoir Parameters: (PhD Dissertation). Qingdao. Ocean University of China, 2008. the

Claims (1)

1. A pre-stack nonlinear inversion method based on particle swarm optimization algorithm, using amplitude-preserved high-resolution seismic data for nonlinear inversion, using updated reflection coefficient formulas in the process of solving nonlinear problems, and using particle swarm optimization The algorithm is optimized, and the common reflection point gather in the reflection angle domain obtained by the amplitude preservation process is used as input. The specific steps are: Step 1: For any location underground, firstly, the value ranges of the compressional wave velocity reflection coefficient, the shear wave velocity reflection coefficient, the density reflection coefficient, and the shear to longitudinal wave velocity ratio are given, the maximum number of iterations and the size of the population are set, and the weighting factor and Learning factor, setting the allowable error range; Step 2: Initialize the longitudinal-wave velocity reflection coefficient, shear-wave velocity reflection coefficient, density reflection coefficient, and shear-to-pillar velocity ratio of each particle with random numbers, and set the initial velocity to zero, so that the individual optimal value and initial value of the particle The same, and calculate the fitness of each particle, select the particle with the smallest fitness as the current global optimal solution; Step 3: Compare whether the fitness of the current global optimum has reached the tolerance range of the error, if it reaches the tolerance range of the error, stop the calculation, output the global optimum as the inversion result, otherwise go to step 4; Step 4: Update the velocity of the particle, and update the position of the particle with the new particle velocity. The number of iterations is increased once. If the number of iterations exceeds the maximum number of iterations, stop the calculation and output the global optimum; Step 5: Recalculate the fitness of the particles, and check whether the fitness of each particle is less than the fitness before updating. If it is less than, update the individual optimum of the particle to the updated position of the particle; Step 6: Compare the new individual optimum, select the individual optimum with the smallest fitness as the global optimum, judge whether the global optimum reaches the error tolerance range, and output the global optimum as the inversion result if the error tolerance range is reached, otherwise return to Step 4, until all calculation points underground are processed; The compressional wave velocity reflection coefficient, the shear wave velocity reflection coefficient, the density reflection coefficient, and the shear to longitudinal wave velocity ratio in each described step are defined in the following manner: In geophysics, the reflection coefficient of a continuous function x(t) is defined as:

Sample a continuous function x(t) into a discrete function with a sampling interval of dt: where x can be P-wave velocity α, S-wave velocity β, and density ρ; subscript 1 indicates the parameters of the underlying formation, and subscript 2 indicates the parameters of the overlying formation; define the ratio of S-wave velocity to P-wave:   

It is characterized in that, the computing method of the particle fitness in described step 2,3,5 is as follows: Define the matrix:

A, B, C, D and K in the formula are defined in the following way:

where i ₁ is the incident angle; Define the matrix:

From the Zoplitz equation, the solution of the following form can be obtained:

From this expression form, the expression of particle fitness can be obtained as:

In the formula, E is particle fitness, w(θ) is seismic wavelet, and d(θ) is seismic observation data.

Images