CN103823239A - Frequency domain optimization mixture staggered grid finite difference forward modeling method - Google Patents

Abstract

A frequency domain optimization mixture staggered grid finite difference forward modeling method comprises the following steps: 1) providing a time-domain two-dimension sound wave equation; 2) eliminating artificial boundary reflection and obtaining a time-domain two-dimension sound wave equation with a perfectly-matched layer boundary condition; 3) performing Fourier transform on time variables at the two sides of the equation and obtaining a frequency-domain sound wave equation; 4) carrying out finite difference discretization on the frequency-domain sound wave equation with the perfectly-matched layer boundary condition based on a conventional staggered grid and obtaining a finite difference discretization format; 5) carrying out finite difference discretization on the frequency-domain sound wave equation with the perfectly-matched layer boundary condition based on a rotation staggered grid and obtaining a finite difference discretization format; 6) performing optimization mixing on the conventional staggered grid and the rotation staggered grid, grid difference item being weighted average of the grid difference item in the two grid systems, and quality acceleration item being weighted average of a center point and eight points around the center point; and 7) under the criterion of the minimum phase velocity error, calculating an optimized coefficient. The weighting coefficient enables frequency dispersion error due to the finite difference discretization to be the smallest, and the precision of the frequency-domain forward modeling is greatly improved.

Description

Frequency field optimization mixes staggered-mesh finite difference analogy method

Technical field

The invention belongs to seismic event forward simulation technology field, relate in particular to a kind of frequency field optimization by the conventional staggered-mesh of frequency field and the two nested grid differential system optimal combinations of frequency field rotationally staggered grid and mix staggered-mesh finite difference analogy method, be intended to improve frequency field forward simulation precision and counting yield.

Background technology

The advantages such as in seismic event forward simulation prior art field, compared with time domain finite-difference forward modeling, it is high that frequency field finite-difference forward modeling has many big guns analog computation efficiency, easy simulated formation absorption effect.Wave equation migration and Full wave shape inverting only need limited several frequencies conventionally, and at this moment frequency field forward simulation more shows the advantage that it is unique.

Analysis frequency territory forward simulation equation can find, equation can be reduced to quality acceleration term and map declination subitem (Laplace operator) sum is zero; In order to improve the precision of frequency field finite-difference forward modeling, conventional staggered-mesh and rotationally staggered grid are introduced to quality acceleration term prioritization scheme, be expressed as central point and its weighted mean of 4 around by quality acceleration term, conventional staggered-mesh and optimize rotationally staggered grid is optimized respectively; Conventional staggered-mesh, rotationally staggered grid, the conventional staggered-mesh of optimization and optimization rotationally staggered grid are 5 difference schemes.

But identical grid difference system, frequency field finite-difference forward modeling is lower than the precision of time domain finite-difference forward modeling, therefore improve precision and become current frequency field forward simulation important research content and be badly in need of solve technical barrier.

Summary of the invention

Object of the present invention is intended to overcome above-mentioned prior art grid difference system than the low defect of time domain finite-difference forward modeling precision, and provide a kind of, conventional staggered-mesh and rotationally staggered grid are optimized to mixing, mixing staggered-mesh finite difference analogy method is optimized in forming frequency territory, and the method can improve high precision and the counting yield of frequency field forward simulation.

Frequency field optimization of the present invention mixes staggered-mesh finite difference analogy method, comprises the steps:

A, provide time domain two dimension ACOUSTIC WAVE EQUATION, single order stress-velocity equation group form;

B, the reflection of elimination Artificial Boundaries, add complete matching layer boundary condition to equation, obtains the time domain two dimension ACOUSTIC WAVE EQUATION with complete matching layer boundary condition;

C, the time domain two dimension ACOUSTIC WAVE EQUATION the right and left time variable with complete matching layer boundary condition is carried out to Fourier transform obtain the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition;

D, that the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition is carried out to finite difference according to conventional staggered-mesh is discrete, obtains conventional staggered-mesh finite difference discrete scheme;

E, the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition is carried out to the discrete rotationally staggered grid finite difference discrete scheme that obtains of finite difference according to rotationally staggered grid;

F, conventional staggered-mesh and rotationally staggered grid two are overlapped to net system carry out hybrid optimization, map declination subitem-Laplace operator is expressed as the weighted mean of map declination subitem-Laplace operator in conventional staggered-mesh and rotationally staggered grid two nested grid systems, quality acceleration term is expressed as central point and its weighted mean of 8 around, is optimized and mixes staggered-mesh finite difference discrete scheme;

The expression formula of G, derivation phase velocity frequency dispersion error, phase velocity frequency dispersion error is the function of weighting coefficient, under the criterion of phase velocity frequency dispersion error minimum, utilizes optimizing algorithm, as method of conjugate gradient, simulated annealing, to ask for optimization coefficient.

Optimization mixing staggered-mesh in described step F is 9 difference schemes.

The present invention compared with prior art remarkable result is:

(1) conventional staggered-mesh and rotationally staggered grid are optimized mixing by the present invention, form to optimize and mix staggered-mesh, obtain 9 difference algorithms, the wherein weighted mean of Laplace operator in two nested grid systems for Laplace operator, quality acceleration term is expressed as central point and its weighted mean of 8 around, and the method has greatly improved the precision of frequency field finite-difference forward modeling;

(2) quality acceleration term weighting coefficient of the present invention optimization can reduce grid anisotropy, improve the precision of forward simulation, optimize conventional staggered-mesh less than conventional staggered-mesh frequency dispersion error, forward simulation precision is high, same optimization rotationally staggered grid is less than rotationally staggered grid frequency dispersion error, and forward simulation precision is high; The optimization of map declination subitem (Laplace operator) weighting coefficient can further reduce frequency dispersion error, optimizes and mixes staggered-mesh anisotropy and frequency dispersion error minimum, and forward simulation precision is the highest;

(3) the present invention optimizes weighting coefficient and makes the discrete frequency dispersion error minimum causing of finite difference, has improved greatly the precision of frequency field forward simulation.

Accompanying drawing explanation

Fig. 1 is the conventional staggered-mesh schematic diagram of the present invention

Fig. 2 is rotationally staggered grid schematic diagram of the present invention

Fig. 3 is optimization mixing staggered-mesh schematic diagram provided by the invention

Fig. 4 a is the conventional staggered-mesh phase velocities dispersion curve of the present invention schematic diagram

Fig. 4 b is the conventional staggered-mesh group velocity dispersion of the present invention curve synoptic diagram

Fig. 5 a is that the present invention optimizes conventional staggered-mesh phase velocities dispersion curve schematic diagram

Fig. 5 b is that the present invention optimizes conventional staggered-mesh group velocity dispersion curve synoptic diagram

Fig. 6 a is rotationally staggered grid phase velocities dispersion curve schematic diagram of the present invention

Fig. 6 b is rotationally staggered grid group velocity dispersion curve synoptic diagram of the present invention

Fig. 7 a is that the present invention optimizes rotationally staggered grid phase velocities dispersion curve schematic diagram

Fig. 7 b is that the present invention optimizes rotationally staggered grid group velocity dispersion curve synoptic diagram

Fig. 8 a is that the present invention optimizes mixing staggered-mesh phase velocities dispersion curve schematic diagram

Fig. 8 b is that the present invention optimizes mixing staggered-mesh group velocity dispersion curve synoptic diagram

Fig. 9 a is finite-difference forward modeling experiment 4Hz real part wave field schematic diagram of the present invention

Fig. 9 b is finite-difference forward modeling experiment 4Hz imaginary part wave field schematic diagram of the present invention

Fig. 9 c is finite-difference forward modeling experiment 10Hz real part wave field schematic diagram of the present invention

Fig. 9 d is finite-difference forward modeling experiment 10Hz imaginary part wave field schematic diagram of the present invention

Figure 10 a is finite-difference forward modeling experiment layered medium model schematic diagram of the present invention

Figure 10 b is finite-difference forward modeling experiment 4Hz real part wave field schematic diagram of the present invention

Figure 10 c is finite-difference forward modeling experiment 4Hz imaginary part wave field schematic diagram of the present invention

Figure 10 d is finite-difference forward modeling experiment 10Hz real part wave field schematic diagram of the present invention

Figure 10 e is finite-difference forward modeling experiment 10Hz imaginary part wave field schematic diagram of the present invention

Above-mentioned Fig. 1-Figure 10 is the accompanying drawing that computing machine painting software suffer draws

Embodiment

Below in conjunction with accompanying drawing and instantiation, the invention will be further described

Frequency field optimization mixes staggered-mesh finite difference analogy method, by acoustic wave equation in time domain (single order stress-velocity equation group), add complete matching layer boundary condition, again system of equations the right and left time variable is carried out to Fourier transform and obtain frequency field ACOUSTIC WAVE EQUATION, according to conventional staggered-mesh (as shown in Figure 1) and rotationally staggered grid (as shown in Figure 2), frequency field ACOUSTIC WAVE EQUATION is carried out to difference discrete respectively, this discretize difference equation can be noted by abridging and be equalled zero for quality acceleration term and map declination subitem sum, quality acceleration term is expressed as to central point and the weighted mean of 4 around, conventional staggered-mesh and optimization rotationally staggered grid are optimized, quality acceleration term is expressed as to central point and its weighted mean of 8 around, map declination subitem (Laplace operator) is expressed as the weighted mean of map declination subitem (Laplace operator) in conventional staggered-mesh and rotationally staggered grid two nested grid systems, is optimized and mixes staggered-mesh (as shown in Figure 3).The method implementation procedure comprises the following steps:

1, provide time domain two dimension ACOUSTIC WAVE EQUATION, single order stress-velocity equation group forms;

$\left\{\begin{array}{c}\frac{\∂P(x,z,t)}{\∂t}=-\mathrm{\κ}(x,z)(\frac{\∂U(x,z,t)}{\∂x}+\frac{\∂V(x,z,t)}{\∂z})\\ \frac{\∂U(x,z,t)}{\∂t}-b(x,z)\frac{\∂P(x,z,t)}{\∂x}+S_x(x,z,t)\\ \frac{\∂V(x,z,t)}{\∂t}=-b(x,z)\frac{\∂P(x,z,t)}{\∂z}+S_z(x,z,t)\end{array}\right.----\left(1\right)$

P in formula (x, z, t) represents hydrostatic force, U (x, z, t) represent the axial speed of x, V (x, z, t) represents the axial speed of z, κ (x, z) represent bulk modulus, the inverse that b (x, z) is density, S _x(x, z, t) and S _z(x, z, t) represents respectively x direction of principal axis and the axial focus item of z.

2, for eliminating Artificial Boundaries reflection, equation is added to complete matching layer boundary condition, obtain the time domain two dimension ACOUSTIC WAVE EQUATION with complete matching layer boundary condition;

$\left\{\begin{array}{c}\frac{\∂P_x(x,z,t)}{\∂t}+d_x(x,z)P_x(x,z,t)=-\mathrm{\κ}(x,z)\frac{\∂U(x,z,t)}{\∂x}\\ \frac{\∂P_z(x,z,t)}{\∂t}+d_z(x,z)P_z(x,z,t)=-\mathrm{\κ}(x,z)\frac{\∂V(x,z,t)}{\∂z}\\ \frac{\∂U(x,z,t)}{\∂t}+d_x(x,z)U(x,z,t)=-b(x,z)(\frac{\∂P_x(x,z,t)}{\∂x}+\frac{\∂P_z(x,z,t)}{\∂x})+S_x(x,z,t)\\ \frac{\∂V(x,z,t)}{\∂t}+d_z(x,z)V(x,z,t)=-b(x,z)(\frac{\∂P_x(x,z,t)}{\∂z}+\frac{\∂P_z(x,z,t)}{\∂z})+S_z(x,z,t)\end{array}\right.---\left(2\right)$

D in formula _x(x, z) is the one dimension function of x coordinate, same d _z(x, z) is also the one dimension function of z coordinate, and their functional form has consistance, conventionally gets

l represents the thickness of matching layer, and x is a local coordinate value, represents the distance of net point to matching layer inward flange.

3, equation (2) the right and left time variable is carried out Fourier transform and is obtained the frequency field ACOUSTIC WAVE EQUATION of complete matching layer (PML) boundary condition;

$\left\{\begin{array}{c}\frac{-\mathrm{i\ω}{\mathrm{\ξ}}_x(x,z)}{\mathrm{\κ}(x,z)}{P}_x(x,z,\mathrm{\ω})=\frac{\∂U(x,z,\mathrm{\ω})}{\∂x}\\ \frac{-\mathrm{i\ω}{\mathrm{\ξ}}_z(x,z)}{\mathrm{\κ}(x,z)}{P}_z(x,z,\mathrm{\ω})=\frac{\∂U(x,z,\mathrm{\ω})}{\∂z}\\ U(x,z,\mathrm{\ω})=-\frac{b(x,z)}{\mathrm{i\ω}{\mathrm{\ξ}}_x(x,z)}\frac{\∂P(x,z,\mathrm{\ω})}{\∂x}+\frac{1}{\mathrm{i\ω}{\mathrm{\ξ}}_x(x,z)}{S}_x(x,z,\mathrm{\ω})\\ V(x,z,\mathrm{\ω})=-\frac{b(x,z)}{\mathrm{i\ω}{\mathrm{\ξ}}_z(x,z)}\frac{\∂P(x,z,\mathrm{\ω})}{\∂z}+\frac{1}{\mathrm{i\ω}{\mathrm{\ξ}}_z(x,z)}{S}_z(x,z,\mathrm{\ω})\end{array}\right.---\left(3\right)$

In formula

p _x(x, z, ω), P _z(x, z, ω), U (x, z, ω), V (x, z, ω), S _x(x, z, ω) and S _z(x, z, ω) is respectively P _x(x, z, t), P _z(x, z, t), U (x, z, t), V (x, z, t), S _x(x, z, t) and S _zthe Fourier transform of (x, z, t).

4, equation (3) is carried out according to conventional staggered-mesh (parameter-definition is shown in Fig. 1) to finite difference is discrete obtains conventional staggered-mesh finite difference discrete scheme, not CONSIDERING BOUNDARY CONDITIONS, discretize difference equation can be abbreviated as:

$\frac{{\mathrm{\ω}}^2}{{\mathrm{\κ}}_{i,j}}{P}_{i,j}+{\mathrm{\Γ}}_{i,j}=S_{i,j}---\left(4\right)$

In formula

be called quality acceleration term, ${\mathrm{\Γ}}_{i,j}=\frac{1}{{\mathrm{\Δ}}^2}(P_{i,j-1}+P_{i-1,j}+P_{i+1,j}+P_{i,j+1}-4P_{i,j})$ Be called map declination subitem (Laplace operator).

5, equation (3) is carried out to the discrete rotationally staggered grid finite difference discrete scheme that obtains of finite difference according to rotationally staggered grid (parameter-definition is shown in Fig. 2), not CONSIDERING BOUNDARY CONDITIONS, equation can be abbreviated as:

$\frac{{\mathrm{\ω}}^2}{{\mathrm{\κ}}_{i,j}}{P}_{i,j}+R_{i,j}=S_{i,j}---\left(5\right)$

In formula

be called quality acceleration term, $R_{\mathrm{ij}}=\frac{1}{2{\mathrm{\Δ}}^2}(P_{i-1,j-1}+P_{i+1,j-1}+P_{i-1,j+1}+P_{i+1,j+1}-4P_{i,j})$ Be called map declination subitem (Laplace operator).

6, conventional staggered-mesh and rotationally staggered grid two are overlapped to net system and carry out hybrid optimization, map declination subitem (Laplace operator) is expressed as the weighted mean of map declination subitem (Laplace operator) in two nested grid systems, and quality acceleration term is expressed as central point and the weighted mean of 8 around it

$\begin{array}{c}\frac{a_1}{{\mathrm{\&Delta;}}^2}(P_{i,j-1}+P_{i-1,j}+P_{i+1,j}+P_{i,j+1}-4P_{i,j})+\frac{1-a_1}{2{\mathrm{\&Delta;}}^2}(P_{i-1,j-1}+P_{i+1,j-1}+P_{i-1,j+1}+P_{i+1,j+1}-4P_{i,j})+\\ \frac{{\mathrm{\&omega;}}^2}{v^2}[c_1{P}_{i,j}+c_2(P_{i.j-1}+P_{i-1,j}+P_{i+1,j}+P_{i+1,j+1})+\frac{1-c_1-4*{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="398519DEST\_PATH\_1403041329201.png,603236DEST\_PATH\_1403041329212.png"\; state="\{\{state\}\}">c</figure-callout>}}_2}{4}(P_{i-1,j-1}+P_{i+1,j-1}+P_{i-1,j+1}+P_{i+1,j+1})]=S_{\mathrm{ij}}\end{array}---\left(6\right)$

The optimization mixing staggered-mesh finite difference method that the present invention that Here it is obtains.

7, can find out conventional staggered-mesh, rotationally staggered grid, optimizes conventional staggered-mesh and optimizes rotationally staggered grid and can regard as to optimize and mix in staggered-mesh weighting coefficient and get a specific class value and obtain.Weighting coefficient is tried to achieve under the criterion of phase velocity frequency dispersion error minimum, and the expression formula of the phase velocity frequency dispersion error that needs here to derive, then can ask for optimization weighting coefficient by optimization algorithm.Weighting coefficient is asked for and be the results are shown in following table.

Table 1 staggered-mesh optimization coefficient

Conventional staggered-mesh and rotationally staggered grid are optimized mixing by the present invention, form to optimize and mix staggered-mesh, obtain 9 difference algorithms, the wherein weighted mean of Laplace operator in two nested grid systems for Laplace operator, quality acceleration term is expressed as central point and its weighted mean of 8 around, and the method has greatly improved the precision of frequency field finite-difference forward modeling.

Phase velocity frequency dispersion error is the function of weighting coefficient, under the criterion of phase velocity frequency dispersion error minimum, utilizes optimization algorithm can ask for one group of optimization weighting coefficient.As shown in Figure 4, Fig. 4 a, 4b have provided phase velocity and the group velocity dispersion curve of conventional staggered-mesh finite difference simulation; As shown in Figure 5, Fig. 5 a, 5b have provided phase velocity and the group velocity dispersion curve of optimizing conventional staggered-mesh finite difference simulation; As shown in Figure 6, Fig. 6 a, 6b have provided phase velocity and the group velocity dispersion curve of rotationally staggered grid finite-difference forward modeling; As shown in Figure 7, Fig. 7 a, 7b have provided phase velocity and the group velocity dispersion curve of optimizing rotationally staggered grid finite-difference forward modeling; As shown in Figure 8, Fig. 8 a, 8b have provided and have optimized the phase velocity and the group velocity dispersion curve that mix staggered-mesh finite difference simulation.Can find out that from optimum results the optimization of quality acceleration term weighting coefficient can reduce grid anisotropy, improve the precision of forward simulation, optimize conventional staggered-mesh less than conventional staggered-mesh frequency dispersion error, forward simulation precision is high, same optimization rotationally staggered grid is less than rotationally staggered grid frequency dispersion error, and forward simulation precision is high; The optimization of map declination subitem (Laplace operator) weighting coefficient can further reduce frequency dispersion error, optimizes and mixes staggered-mesh grid anisotropy and frequency dispersion error minimum, and forward simulation precision is the highest.

What the present invention provided is based on ACOUSTIC WAVE EQUATION derivation, can the method be generalized to equations for elastic waves according to above-mentioned thinking.

Application example

In order to prove to optimize the validity and the essential characteristic of showing frequency field wave field of mixing staggered-mesh finite difference analogy method, the invention provides two examples, one is uniform dielectric forward modeling simulative example, and another is layered medium forward modeling simulative example.Uniform dielectric forward simulation is referring to Fig. 9, Fig. 9 a, 9b have provided respectively 4Hz wave field real part and imaginary part, Fig. 9 c, 9d have provided respectively 10Hz wave field real part and imaginary part, can find out that the single-frequency wave field that uniform dielectric medium frequency territory forward simulation obtains is the concentric circles around focus.Referring to Figure 10, provide layered medium forward modeling analog result, Figure 10 a has provided layered medium illustraton of model, Figure 10 b, 10c have provided respectively 4Hz wave field real part and imaginary part, Figure 10 d, 10e have provided respectively 10Hz wave field real part and imaginary part, can find out that the single-frequency wave field that layered medium medium frequency territory forward simulation obtains is no longer complete concentric circles.When can Inference Model more complicated, frequency field single-frequency wave field will become more complicated, be not easy to find out the precision of forward simulation from frequency field single-frequency wave field, the precision that contrasts in detail various difference method forward simulations can contrast the dispersion curve of each grid difference system.

Claims (2)

1. frequency field optimization mixes staggered-mesh finite difference analogy method, it is characterized in that: the method comprises the steps:

A, provide time domain two dimension ACOUSTIC WAVE EQUATION, single order stress-velocity equation group form;

B, the reflection of elimination Artificial Boundaries, add complete matching layer boundary condition to equation, obtains the time domain two dimension ACOUSTIC WAVE EQUATION with complete matching layer boundary condition;

C, the time domain two dimension ACOUSTIC WAVE EQUATION the right and left time variable with complete matching layer boundary condition is carried out to Fourier transform obtain the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition;

D, that the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition is carried out to finite difference according to conventional staggered-mesh is discrete, obtains conventional staggered-mesh finite difference discrete scheme;

E, the frequency field ACOUSTIC WAVE EQUATION with complete matching layer boundary condition is carried out to the discrete rotationally staggered grid finite difference discrete scheme that obtains of finite difference according to rotationally staggered grid;

F, conventional staggered-mesh and rotationally staggered grid two are overlapped to net system carry out hybrid optimization, map declination subitem-Laplace operator is expressed as the weighted mean of map declination subitem-Laplace operator in conventional staggered-mesh and rotationally staggered grid two nested grid systems, quality acceleration term is expressed as central point and its weighted mean of 8 around, is optimized and mixes staggered-mesh finite difference discrete scheme;

The expression formula of G, derivation phase velocity frequency dispersion error, phase velocity frequency dispersion error is the function of weighting coefficient, under the criterion of phase velocity frequency dispersion error minimum, utilizes optimizing algorithm, as method of conjugate gradient, simulated annealing, to ask for optimization coefficient.

2. frequency field optimization according to claim 1 mixes staggered-mesh finite difference analogy method, it is characterized in that, the optimization mixing staggered-mesh in described step F is 9 difference schemes.

Images