CN103698814A - Implementation method for mixed absorbing boundary condition applied to variable density acoustic wave equation - Google Patents

Abstract

The invention relates to an implementation method for a mixed absorbing boundary condition applied to a variable density acoustic wave equation. The implementation method comprises the steps of dividing a whole simulation region into an internal region and a boundary region, calculating the wave field value of a sound wave of each absorbing boundary layer by using the variable density acoustic wave equation and a variable density AWWE (Arbitrary Wide-angle wave equation) respectively, and performing linear weighting on the two calculation results to obtain the wave field value of the sound wave of each absorbing boundary layer. On the basis of realizing a single variable density AWWE absorbing boundary layer and in 2D (2 dimensional) and 3D (3 dimensional) space coordinate systems, the wave field value of the sound wave of each absorbing boundary layer is calculated by respectively adopting a combination mode of processing edges of the variable density AWWE and processing angles of a modified 15-degree one-way wave equation, and a combination mode of processing faces of the variable density AWWE, the processing edges of the modified 15-degree one-way wave equation and the processing angles of a 5-degree one-way wave equation, and the mixed absorbing boundary condition of the variable density AWWE is realized. The implementation method can be widely applied to the numerical simulation of constant density and variable density acoustic wave equations.

Description

A kind of implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION

Technical field

The present invention relates to a kind of implementation method of mixed absorbing boundary, particularly about a kind of implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION.

Background technology

In field of seismic exploration, Wave Equation Numerical technology is to understand the important means of seismic wave propagation rule, and it not only can instruct the explanation of seismic data and the design of recording geometry, or the basis of current popular treatment technology reverse-time migration and Full wave shape inverting.Wave Equation Numerical normally carries out in limited region, and the outside of simulated domain exists manual intercept border.If do not carry out special processing, when propagating into cutoff boundary, ripple can produce edge reflection, thus normal reflective information intrusively.Effective disposal route is an energy for traveling-wave field outside the surrounding of simulated domain is utilized absorbing boundary absorption, thus compacting edge reflection.

Current absorbing boundary condition of a great variety, the absorbing boundary based on one way wave equation is the large class of one in absorbing boundary.The thinking of one way ripple absorbing boundary is to utilize one way wave equation to calculate wave field at boundary, makes seismic event have single directivity at boundary, thereby avoids wave field to be reflected back toward internal simulation region.Yet one way ripple absorbing boundary exists the shortcoming of angle limits, early stage 5 ° and 15 ° of low order one way wave equations, can not absorb the energy of wide-angle incident wave and angle point incident wave well.In order to expand one way wave propagation angle, some high-order one way wave equations are developed successively, and wherein that more famous is AWWE(Arbitrary Wide-angle Wave Equation, arbitrarily wide-angle wave equation).AWWE has two obvious advantages, and the one, support high dip angle, its propagation angle approaches 90 °; The 2nd, form is simple, only comprises second-order partial differential coefficient item, and numerical evaluation is convenient.Because these advantages, are successfully used for migration imaging by AWWE, also it is used as to the absorbing boundary in second order Chang Midu ACOUSTIC WAVE EQUATION numerical simulation simultaneously.What in second order Chang Midu ACOUSTIC WAVE EQUATION numerical simulation, utilize is the AWWE absorbing boundary of Chang Midu, and when utilizing the variable density ACOUSTIC WAVE EQUATION of more realistic one-order velocity-pressure form to carry out numerical simulation, the AWWE absorbing boundary of Chang Midu is no longer applicable.

Although AWWE compares with other low order one way wave equations, the effect of its compacting edge reflection has had obvious improvement, but still can not meet the requirement of high precision wave field numerical.A kind of thinking of improving one way ripple absorbing boundary assimilation effect is to mix round trip ripple and one way ripple, and its strategy is multilayer absorption edge interlayer to be set in frontier district form transitional zone, rather than conventional individual layer, and this class absorbing boundary is named as blended absorbent border.Existing mixed absorbing boundary utilizes low order one way wave equation, and its effect that absorbs wide-angle incident wave and angle point incident wave is still not ideal enough, the hypothesis of existing blended absorbent border based on Chang Midu simultaneously, and be not suitable for the ACOUSTIC WAVE EQUATION of variable density.

Summary of the invention

For the problems referred to above, the implementation method that the object of this invention is to provide a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION, the blended absorbent border of realizing by the method can be applicable to the numerical simulation of 2D/3D variable density ACOUSTIC WAVE EQUATION, and its compacting edge reflection is satisfactory for result.

For achieving the above object, the present invention takes following technical scheme: a kind of implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION, it comprises the following steps: 1), when the acoustic wavefield of simulation variable density, whole 2D computer memory is divided into 2D interior zone and comprises single absorption edge interlayer I ₁j ₁k ₁l ₁2D absorbing boundary region; Utilize the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on the outer boundary ABCD of 2D interior zone, according to the acoustic wavefield value on the outer boundary ABCD calculating, utilize the AWWE of the variable density under 2D space coordinates to process the array mode of 15 ° of one way wave equation processing angles amended under rib and 2D space coordinates, the single absorption edge interlayer I in calculating absorbing boundary region ₁j ₁k ₁l ₁on acoustic wavefield value; 2) the single absorption edge interlayer I in 2D absorbing boundary region ₁j ₁k ₁l ₁arranged outside n-1 layer AWWE absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density of the variable density under 2D space coordinates, the mode by multi-ply linear weighting calculates clathrum ω ₂=1,2 ..., the acoustic wavefield value on n; 3) 2D computer memory is expanded to 3D computer memory, on the up, down, left, right, before and after six direction in 3D internal calculation region, single absorbing boundary is set simultaneously, according to the acoustic wavefield value on the outer boundary ABCD of 3D interior zone, utilize amended 15 ° of one way wave equations under the AWWE treated side, 2D space coordinates of the variable density under 3d space coordinate system to process the array mode of 5 ° of one way wave equation processing angles under rib and 2D space coordinates, calculate the acoustic wavefield value on the single absorption edge interlayer of 3D interior zone six direction; 4) outside of the single absorption edge interlayer simultaneously arranging on the up, down, left, right, before and after six direction of 3D interior zone all arranges multilayer absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density of the variable density under 3d space coordinate system, the mode by multi-ply linear weighting calculates clathrum ω ₃=1,2 ..., the acoustic wavefield value on n, the AWWE blended absorbent border of realizing 3D variable density.

In described step 1), utilize the AWWE and amended 15 ° of modes that one way wave equation is combined of the variable density under 2D space coordinates, calculate the single absorption edge interlayer I in absorbing boundary region ₁j ₁k ₁l ₁on acoustic wavefield value, it comprises the following steps: the outer boundary ABCD of (1) 2D interior zone comprises rib BC, rib CD, rib DA and rib AB, utilizes the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on each rib in outer boundary ABCD; (2) according to the acoustic wavefield value on the rib BC of the outer boundary ABCD of 2D interior zone, rib CD, rib DA and rib AB, utilize respectively the AWWE of descending, right lateral under 2D space coordinates, up and left lateral variable density, calculate lower boundary J ₁k ₁, right margin K ₁l ₁, coboundary L ₁i ₁with left margin I ₁j ₁in four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value; (3) utilize amended 15 ° of one way wave equations under 2D space coordinates, respectively according to additional clathrum E ₁f ₁g ₁h ₁the acoustic wavefield value of mid point M, N, O, P, Q, R, S, T, calculates additional clathrum E ₁f ₁g ₁h ₁in the outer end points W of four rib MN, OP, QR, ST ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value; (4) the single absorption edge interlayer I obtaining according to step (2) ₁j ₁k ₁l ₁in four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value and the end points W that obtains of step (3) ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value, realize the AWWE absorbing boundary condition of the simple layer variable density under 2D space coordinates.

In described step (2), four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on the computing method of acoustic wavefield value specifically comprise the following steps: 1. according to the acoustic wavefield value on rib BC, rib CD, rib DA and rib AB in the outer boundary ABCD of 2D interior zone, the AWWE that utilizes descending, right lateral under 2D space coordinates, up and left lateral variable density, calculates additional clathrum E by finite difference ₁f ₁g ₁h ₁acoustic wavefield value on middle rib MN, rib OP, rib QR, rib ST; 2. according to following formula

$p_{k+1}=2p_{k+\frac{1}{2}}-p_k,$

Respectively by the acoustic wavefield value substitution p on rib BC, rib CD, rib DA and rib AB _k, respectively by the acoustic wavefield value on rib MN, rib OP, rib QR and rib ST the p that correspondence calculates _k+1value be single absorption edge interlayer I ₁j ₁k ₁l ₁middle rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value.

The equation form of the AWWE of the descending and up variable density under described 2D space coordinates is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}\&PlusMinus;{\mathrm{<figure-callout\; id="1"\; label="cH"\; filenames="HDA0000451812910000042.png"\; state="\{\{state\}\}">cH</figure-callout>}}_1\frac{\&PartialD;}{\&PartialD;z})p={\mathrm{\&rho;c}}^2{\mathrm{<figure-callout\; id="2"\; label="H"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_2\frac{\&PartialD;v_x}{\&PartialD;x}\\ \frac{\&PartialD;v_x}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;x}\end{array}\right.,$

In formula, in the equation of the AWWE of descending variable density

get "+", in the equation of the AWWE of up variable density

get "-"; X and z represent the coordinate under 2D space coordinates, and t represents time coordinate; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing;

v _xrepresent the speed component of vibration in the x-direction,

for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; H ₁and H ₂equal representing matrix; X in the equation of the AWWE of descending variable density and z are exchanged, obtain the equation of the AWWE of right lateral variable density; X in the equation of the AWWE of up variable density and z are exchanged, obtain the equation of the AWWE of left lateral variable density.

Under described 2D space coordinates, the form of amended left lateral and 15 ° of one way wave equations of right lateral is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;x\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;\mathrm{<figure-callout\; id="2"\; label="t"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">t</figure-callout>}}^2}+\frac{c}{2}(\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">c</figure-callout>}}^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;x}^2})=0;$

In formula, in 15 ° of one way wave equations of amended left lateral

get "+", in 15 ° of one way wave equations of amended right lateral

get "-"; Change the x in 15 ° of one way wave equations of amended left lateral into z, obtain amended up 15 ° of one way wave equations under 2D space coordinates; Change the x in 15 ° of one way wave equations of amended right lateral into z, obtain amended descending 15 ° of one way wave equations under 2D space coordinates.

Described step 2) in, computing grid layer ω ₂=1,2 ..., the acoustic wavefield value on n, it specifically comprises the following steps: (1) utilizes the ACOUSTIC WAVE EQUATION computing grid layer ω of the variable density under 2D space coordinates ₂=1,2 ..., the acoustic wavefield value on n-1, and by clathrum ω ₂acoustic wavefield value on=n is preset as zero, clathrum ω ₂=1,2 ..., n-1, the acoustic wavefield value on n is designated as

(2) according to the acoustic wavefield value on the outer boundary ABCD of 2D interior zone, utilize the variable density under 2D space coordinates AWWE and

extrapolation calculates clathrum ω ₂acoustic wavefield value v on=1 ₁; Respectively according to the clathrum ω being obtained by step 1) ₂=1,2 ..., the acoustic wavefield value on n-1 utilize the variable density under 2D space coordinates AWWE and

extrapolation calculates clathrum ω ₂=2 ..., n-1, the acoustic wavefield value on n (3) the clathrum ω that will be obtained by step (1) ₂=1,2 ..., n-1, the acoustic wavefield value on n

with the clathrum ω being obtained by step (2) ₂=1,2 ..., n-1, the acoustic wavefield value on n

carry out linear weighted function, obtain clathrum ω after weighting ₂=1,2 ..., n-1, the acoustic wavefield value on n for:

$p_{{\mathrm{\&omega;}}_2}=\frac{{\mathrm{\&omega;}}_2}{n}{v}_{{\mathrm{\&omega;}}_2}+\frac{n-{\mathrm{\&omega;}}_2}{n}{u}_{{\mathrm{\&omega;}}_2},({\mathrm{\&omega;}}_2=\mathrm{1,2},...,n).$

In described step 3), the implementation method of the single absorption edge interlayer on the up, down, left, right, before and after six direction in 3D internal calculation region is identical, wherein calculates the single absorption edge interlayer I of 3D interior zone down direction ₃j ₃k ₃l ₃on acoustic wavefield value comprise the following steps: (1), according to the acoustic wavefield value on the outer boundary ABCD in 3D internal calculation region, utilizes the AWWE of the descending variable density under 3d space coordinate system, calculates additional clathrum

upper region M ₃n ₃o ₃p ₃four interior rib M ₃n ₃, interior rib N ₃o ₃, interior rib O ₃p ₃with interior rib P ₃m ₃on acoustic wavefield value; (2) according to four interior rib M ₃n ₃, interior rib N ₃o ₃, interior rib O ₃p ₃with interior rib P ₃m ₃on acoustic wavefield value, utilize respectively amendedly under 2D space coordinates move ahead, 15 ° of one way wave equations of rear row, left lateral and right lateral, calculate additional clathrum

upper four outer rib V ₃w ₃, outer rib R ₃s ₃, outer rib X ₃z ₃with outer rib T ₃u ₃on acoustic wavefield value; (3) utilize 5 ° of one way wave equations under 2D space coordinates, respectively according to four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value on two adjoint points separately, calculates additional clathrum

upper four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value; (4) according to formula

respectively by the acoustic wavefield value substitution p on each rib in the outer boundary ABCD of 3D interior zone 11 _k, will add clathrum respectively

in acoustic wavefield value substitution on each rib

the p that correspondence calculates _k+1value be the single absorption edge interlayer I of 3D interior zone down direction ₃j ₃k ₃l ₃in acoustic wavefield value on each rib.

In described step (1), the finite difference computation scheme of the AWWE of the descending variable density under 3d space coordinate system is:

$\left\{\begin{array}{c}{p}_{i,j,k+\frac{1}{2}}^{n+1}+{\mathrm{\&alpha;}}_z{p}_{i,j,k+\frac{1}{2}}^{n+1}a=p_{i,j,k+\frac{1}{2}}^n+{\mathrm{\&alpha;}}_{z}a(p_{i,j,k}^{n+1}+p_{i,j,k}^n-p_{i,j,k+\frac{1}{2}}^n)\\ +\mathrm{\&rho;}{\mathrm{cH}}_2[{\mathrm{\&alpha;}}_x({\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_x\right)}_{i-\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}})+{\mathrm{\&alpha;}}_y({\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_y\right)}_{i,j-\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}})]\\ {\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},j,k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i+1,j,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\\ {\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},j,k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i,j+1,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\end{array}\right.,$

In formula, x, y and z all represent the coordinate under 3d space coordinate system; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing;

v _xrepresent the speed component of vibration in the x-direction,

for calculating the auxiliary variable of introducing, v _yrepresent the speed component of vibration in the y-direction,

for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; H ₂for matrix, i, j and k are respectively the discrete grid block node under 3d space coordinate system, the discrete grid block node that n is the time; Matrix a is: a=(L ₁₁, L ₂₁, L ₃₁..., L _m1) ^t; α _x=c Δ t/ Δ x, α _y=c Δ t/ Δ y, α _z=c Δ t/ Δ z, Δ x, Δ y and Δ z be under 3d space coordinate system mesh spacing, the mesh spacing that Δ t is the time.

In described step (2), under 2D space coordinates, amended moving ahead with the form of 15 ° of one way wave equations of rear row is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;y\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;\mathrm{<figure-callout\; id="2"\; label="t"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">t</figure-callout>}}^2}+\frac{c}{2}(\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">c</figure-callout>}}^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;y}^2})=0,$

In formula, amended moving ahead in 15 ° of one way wave equations

get "+", in 15 ° of one way wave equations of amended rear row

get "-".

In described step (3), 5 ° of one way wave equation forms under 2D space coordinates are:

$\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;x}\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;y}+\frac{\sqrt{2}}{c}\frac{\&PartialD;p}{\&PartialD;t}=0.$

The present invention is owing to taking above technical scheme, it has the following advantages: 1, the present invention is owing to utilizing the AWWE of ACOUSTIC WAVE EQUATION and the variable density of variable density, mode by multi-ply linear weighting calculates the acoustic wavefield value on each clathrum of borderline region, therefore the present invention can adapt to the spatial variations of density, thus the effect of suppressing edge reflection in the situation of assurance density with spatial variations.2,, on the basis of the present invention due to the AWWE absorption edge interlayer in single variable density, in all directions of interior zone, multilayer absorption edge interlayer is set simultaneously, so the present invention can further improve the effect of single absorbing boundary lamination edge reflection processed.3, the present invention is due under 2D space coordinates, utilize the AWWE processing rib of variable density and the array mode of amended 15 ° of one way wave equation processing angles, calculate the borderline acoustic wavefield value of blended absorbent, so the present invention can suppress respectively the acoustic reflection from rib and angle effectively.4, the present invention is due under 3d space coordinate system, utilize amended 15 ° of one way wave equations under the AWWE treated side, 2D space coordinates of variable density to process the array mode of 5 ° of one way wave equation processing angles under rib, 2D space coordinates, calculate the borderline acoustic wavefield value of blended absorbent, so the present invention can suppress respectively the acoustic reflection from face, rib and angle effectively.Based on above advantage, the present invention can be widely used in the numerical simulation of Chang Midu and variable density ACOUSTIC WAVE EQUATION.

Accompanying drawing explanation

Fig. 1 is the AWWE absorbing boundary of individual layer variable density in 2D variable density ACOUSTIC WAVE EQUATION numerical simulation of the present invention and the schematic diagram of Corner Treatment thereof; Wherein, border ABCD is the outer boundary of 2D interior zone, and its place clathrum is labeled as ω ₁=0; Border I ₁j ₁k ₁l ₁for the single absorption edge interlayer in 2D absorbing boundary region, its place clathrum is designated as ω ₁=1; Dotted line E ₁f ₁g ₁h ₁the additional clathrum of introducing when utilizing the AWWE of variable density to calculate, is positioned at outer boundary ABCD and the single absorption edge interlayer I of 2D interior zone ₁j ₁k ₁l ₁centre, its place clathrum is labeled as

Fig. 2 is the schematic diagram on the AWWE blended absorbent border of variable density in 2D variable density ACOUSTIC WAVE EQUATION numerical simulation of the present invention; Wherein, the single absorption edge interlayer I based on shown in Fig. 1 ₁j ₁k ₁l ₁, the n layer AWWE absorption edge interlayer arranging in 2D absorbing boundary region, its corresponding clathrum numbering is from inside to outside labeled as ω successively ₂=1,2 ..., n, the clathrum that is about to the outer boundary ABCD place of 2D interior zone is labeled as ω ₂=0, absorption edge interlayer I ₁j ₁k ₁l ₁the clathrum at place is labeled as ω ₂=1, the like, the outermost layer absorption edge interlayer I in 2D absorbing boundary region ₂j ₂k ₂l ₂the clathrum at place is labeled as ω ₂=n;

Fig. 3 is the schematic diagram of the AWWE absorbing boundary of individual layer variable density in 3D variable density ACOUSTIC WAVE EQUATION numerical simulation of the present invention and rib thereof, Corner Treatment; Wherein, ABCD is the lower boundary of 3D interior zone, and the clathrum at its place is labeled as ω ₃=0; I ₃j ₃k ₃l ₃that the clathrum at its place is labeled as ω at the single absorption edge interlayer of the down direction of 3D interior zone ₃=1; E ₃f ₃g ₃h ₃be the additional clathrum of introducing while utilizing the AWWE of 3D variable density to calculate, it is positioned at lower boundary ABCD and the single absorption edge interlayer I of 3D interior zone ₃j ₃k ₃l ₃centre, its place clathrum is labeled as

Fig. 4 is that individual layer AWWE of the present invention and amended 15 ° of these combined strategies of one way wave equation are for the assimilation effect schematic diagram of 2D uniform dielectric acoustic wavefield numerical simulation, in Fig. 4 (a), utilize individual layer AWWE to process rib, amended 15 ° of one way wave equation processing angles; In Fig. 4 (b), utilize the amended 15 ° of one way wave equations of individual layer to process rib, 5 ° of one way wave equation processing angles, white arrow indication is edge reflection;

Fig. 5 be the AWWE blended absorbent border of 2D variable density of the present invention for the assimilation effect schematic diagram of uniform dielectric acoustic wavefield numerical simulation, Fig. 5 (a) is the acoustic wavefield section of 810ms; Fig. 5 (b) is the wave field section of 1.8s, utilizes individual layer AWWE absorbing boundary; Fig. 5 (c) is the wave field section of 1.8s, utilizes 10 layers of AWWE to form blended absorbent border;

Fig. 6 is the inhomogeneous Marmousi rate pattern schematic diagram of 2D that the present invention utilizes, and corresponding inhomogeneous density model is generated by an experimental formula;

Fig. 7 be 2D variable density AWWE blended absorbent of the present invention border for the assimilation effect schematic diagram of Marmousi model acoustic wavefield numerical simulation, wherein focus is near left margin, degree of depth 20m; Fig. 7 (a) is the wave field section of 825ms, without absorbing boundary; Fig. 7 (b) is the wave field section of 825ms, utilizes the AWWE of 10 layers of variable density to form blended absorbent border; Fig. 7 (c) is the wave field section of 1.87s, without absorbing boundary; Fig. 7 (d) is the wave field section of 1.87s, utilizes the AWWE of 10 layers of variable density to form blended absorbent border;

Fig. 8 is the present invention's seismologic record that above focus, 10m place receives in Fig. 6, in Fig. 8 (a) without absorbing boundary; In Fig. 8 (b), utilize the AWWE of 10 layers of variable density to form blended absorbent border;

Fig. 9 be the AWWE blended absorbent border of 3D variable density of the present invention for the assimilation effect schematic diagram of uniform dielectric acoustic wavefield numerical simulation, source depth 20m, is positioned at upper surface center; Fig. 9 (a) is the wave field section of 112.5ms; Fig. 9 (b) is the wave field section of 250ms, utilizes individual layer AWWE absorbing boundary treated side, and amended 15 ° of one way wave equations are processed rib, 5 ° of one way wave equation processing angles, and white arrow indication is edge reflection; Fig. 9 (c) is the wave field section of 250ms, utilizes 10 layers of AWWE to form blended absorbent border; Fig. 9 (d) is the wave field section of 325ms, utilizes individual layer AWWE absorbing boundary treated side, and amended 15 ° of one way wave equations are processed rib, 5 ° of one way wave equation processing angles, and white arrow indication is edge reflection; Fig. 9 (e) is the wave field section of 325ms, utilizes 10 layers of AWWE to form blended absorbent border;

Figure 10 is the schematic diagram that the present invention utilizes the non-homogeneous rate pattern of 3D, and corresponding density model heterogeneous is generated by an experimental formula;

Figure 11 be the AWWE blended absorbent border of 3D variable density of the present invention for the assimilation effect schematic diagram of non-homogeneous rate pattern acoustic wavefield numerical simulation, focus is positioned at upper surface center, degree of depth 20m; Figure 11 (a) is the section of 155ms wave field, without absorbing boundary; Figure 11 (b) is the section of 155ms wave field, utilizes the AWWE blended absorbent border of 10 layers of variable density; Figure 11 (c) is the section of 350ms wave field, and without absorbing boundary, white arrow indication is edge reflection; Figure 11 (d) is the section of 350ms wave field, utilizes the AWWE of 10 layers of variable density to form blended absorbent border.

Embodiment

Below in conjunction with drawings and Examples, the present invention is described in detail.

The implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION of the present invention comprises the following steps:

1) as shown in Figure 1, when the acoustic wavefield of simulation variable density, by whole 2D(two dimension) computer memory is divided into 2D interior zone 1 and comprises single absorption edge interlayer I ₁j ₁k ₁l ₁2D absorbing boundary region 2; Utilize the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on the outer boundary ABCD of 2D interior zone 1; According to the acoustic wavefield value on the outer boundary ABCD of the 2D interior zone 1 calculating, utilize the AWWE of the variable density under 2D space coordinates to process the array mode of amended 15 ° of one way wave equation processing angles under rib and 2D space coordinates, the single absorption edge interlayer I in calculating absorbing boundary region 2 ₁j ₁k ₁l ₁on acoustic wavefield value, it specifically comprises the following steps:

(1) the outer boundary ABCD of 2D interior zone 1 comprises rib BC, rib CD, rib DA and rib AB, utilizes the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on each rib in outer boundary ABCD.

(2) single absorption edge interlayer I ₁j ₁k ₁l ₁comprise lower boundary J ₁k ₁, right margin K ₁l ₁, coboundary L ₁i ₁with left margin I ₁j ₁, according to the acoustic wavefield value on the rib BC of the outer boundary ABCD of 2D interior zone 1, rib CD, rib DA and rib AB, utilize respectively the AWWE of descending, right lateral under 2D space coordinates, up and left lateral variable density, calculate lower boundary J ₁k ₁, right margin K ₁l ₁, coboundary L ₁i ₁with left margin I ₁j ₁in four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value.

Due to four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on the account form of acoustic wavefield value identical, following rib M only ₁n ₁on the example that is calculated as of acoustic wavefield value describe, it specifically comprises the following steps:

1. according to the acoustic wavefield value on rib BC in the outer boundary ABCD of 2D interior zone 1, utilize the AWWE of the descending variable density under 2D space coordinates, by finite difference, calculate additional clathrum E ₁f ₁g ₁h ₁acoustic wavefield value on middle rib MN.

The equation form of the AWWE of the descending variable density under 2D space coordinates is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}+{\mathrm{<figure-callout\; id="1"\; label="cH"\; filenames="HDA0000451812910000042.png"\; state="\{\{state\}\}">cH</figure-callout>}}_1\frac{\&PartialD;}{\&PartialD;z})p={\mathrm{\&rho;c}}^2{\mathrm{<figure-callout\; id="2"\; label="H"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_2\frac{\&PartialD;v_x}{\&PartialD;x}\\ \frac{\&PartialD;v_x}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;x}\end{array}\right.---\left(1\right)$

In formula (1), x and z represent the coordinate under 2D space coordinates, and t represents time coordinate; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing;

v _xrepresent the speed component of vibration in the x-direction, for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; Matrix H ₁and H ₂be respectively:

$\left\{\begin{array}{c}{H}_1=\mathrm{<figure-callout\; id="2"\; label="LD\; H"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">LD</figure-callout>}& \\ H_{}{}_2={\mathrm{L\&Lambda;}}_2\end{array}\right.---\left(2\right)$

$\left\{\begin{array}{c}L={({\mathrm{\&Lambda;}}_1+{\mathrm{\&Lambda;}}_2)}^{-1}\\ D=d^{T}d\end{array}\right.---\left(3\right)$

In formula (2), d=(1,0 ..., 0) _{1 * m}; In formula (2) and formula (3), matrix Λ ₁and Λ ₂by reference velocity c _i(i=1,2 ..., m) and the matrix that forms of background velocity c, it is respectively:

The finite difference computation scheme of the AWWE of the descending variable density under 2D space coordinates is:

$\left\{\begin{array}{c}{p}_{i,k+\frac{1}{2}}^{n+1}+{\mathrm{\&alpha;}}_z{p}_{i,k+\frac{1}{2}}^{n+1}a=p_{i,k+\frac{1}{2}}^n+{\mathrm{\&alpha;}}_z(p_{i,k}^{n+1}+p_{i,k}^n-p_{i,k+\frac{1}{2}}^n)a+{\mathrm{\&rho;}}_{i,k+\frac{1}{2}}{c}_{i,k+\frac{1}{2}}{H}_2\left[{\mathrm{\&alpha;}}_x({\left(v_x\right)}_{i+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_x\right)}_{i-\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}})\right]\\ {\left(v_x\right)}_{i+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_x\right)}_{i+\frac{1}{2},k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i+1,k+\frac{1}{2}}^n-p_{i,k+\frac{1}{2}}^n)\end{array}\right.---\left(6\right)$

In formula (6), i and k are the discrete grid block node under 2D space coordinates, the discrete grid block node that n is the time, and matrix a is: a=(L ₁₁, L ₂₁, L ₃₁..., L _m1) ^t; α _x=c Δ t/ Δ x, α _z=c Δ t/ Δ z, Δ x and Δ z are respectively the mesh spacing of x direction and z direction under 2D space coordinates, the mesh spacing that Δ t is the time.

2. according to following formula

$p_{k+1}=2p_{k+\frac{1}{2}}-p_k---\left(7\right)$

By the acoustic wavefield value substitution p on rib BC _k, by the acoustic wavefield value substitution on rib MN

calculate p _k+1value, p _k+1value be single absorption edge interlayer I ₁j ₁k ₁l ₁middle rib M ₁n ₁on acoustic wavefield value.

Utilize and calculate lower rib M ₁n ₁on the identical processing mode of acoustic wavefield value, respectively according to the acoustic wavefield value on rib CD, the DA of the outer boundary ABCD of 2D interior zone 1 and AB, the AWWE equation that utilizes right lateral under 2D space coordinates, up, left lateral variable density, calculates additional clathrum E ₁f ₁g ₁h ₁rib OP, QR and the acoustic wavefield value on ST, utilize formula (7) further to obtain single absorption edge interlayer I ₁j ₁k ₁l ₁in other three rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value.

The equation form of the AWWE of the right lateral under 2D space coordinates and left lateral variable density is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}\&PlusMinus;{\mathrm{<figure-callout\; id="1"\; label="cH"\; filenames="HDA0000451812910000042.png"\; state="\{\{state\}\}">cH</figure-callout>}}_1\frac{\&PartialD;}{\&PartialD;x})p={\mathrm{\&rho;c}}^2{\mathrm{<figure-callout\; id="2"\; label="H"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_2\frac{\&PartialD;v_z}{\&PartialD;z}\\ \frac{\&PartialD;v_z}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;z}\end{array}\right.---\left(8\right)$

In formula (8), in the equation of the AWWE of right lateral variable density

get "+", in the equation of the AWWE of left lateral variable density

get "-".

The equation form of the AWWE of the up variable density under 2D space coordinates is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}-{\mathrm{<figure-callout\; id="1"\; label="cH"\; filenames="HDA0000451812910000042.png"\; state="\{\{state\}\}">cH</figure-callout>}}_1\frac{\&PartialD;}{\&PartialD;z})p={\mathrm{\&rho;c}}^2{\mathrm{<figure-callout\; id="2"\; label="H"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_2\frac{\&PartialD;v_x}{\&PartialD;x}\\ \frac{\&PartialD;v_x}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;x}\end{array}\right.---\left(9\right)$

(3) be effectively to suppress the reflection of angle point, utilize amended 15 ° of one way wave equations under 2D space coordinates, respectively according to additional clathrum E ₁f ₁g ₁h ₁the acoustic wavefield value of mid point M, N, O, P, Q, R, S, T, calculates additional clathrum E ₁f ₁g ₁h ₁in the outer end points W of four rib MN, rib OP, rib QR and rib ST ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value, it specifically comprises:

1. utilize the 15 ° of one way wave equations of left lateral under amended 2D space coordinates, the acoustic wavefield value of ordering according to M point and R respectively, calculates left end point W ₁and W ₆acoustic wavefield value; Utilize the 15 ° of one way wave equations of right lateral under amended 2D space coordinates, the acoustic wavefield value of ordering according to N point and Q respectively, calculates right endpoint W ₂and W ₅acoustic wavefield value.

The form of the left lateral under amended 2D space coordinates and 15 ° of one way wave equations of right lateral is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;x\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;\mathrm{<figure-callout\; id="2"\; label="t"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">t</figure-callout>}}^2}+\frac{c}{2}(\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">c</figure-callout>}}^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;x}^2})=0---\left(10\right)$

In formula (10), in 15 ° of one way wave equations of amended left lateral

get "+", in 15 ° of one way wave equations of amended right lateral

get "-".

2. 1. similar with step, utilize the up 15 ° of one way wave equations under amended 2D space coordinates, the acoustic wavefield value of ordering according to S point and P respectively, extrapolation calculates upper extreme point W ₇and W ₄acoustic wavefield value; Utilize the descending 15 ° of one way wave equations under amended 2D space coordinates, the acoustic wavefield value of ordering according to T point and O respectively, extrapolation calculates lower extreme point W ₈and W ₃acoustic wavefield value.

The form of the 15 ° of one way wave equations of uplink and downlink under amended 2D space coordinates is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;z\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;\mathrm{<figure-callout\; id="2"\; label="t"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">t</figure-callout>}}^2}+\frac{c}{2}(\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">c</figure-callout>}}^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;z}^2})=0---\left(11\right)$

In formula (11), in amended up 15 ° of one way wave equations

get "+", in amended descending 15 ° of one way wave equations get "-".

(4) the single absorption edge interlayer I obtaining according to step (2) ₁j ₁k ₁l ₁in four rib M ₁n ₁, O ₁p ₁, Q ₁r ₁, S ₁t ₁on acoustic wavefield value and the end points W that obtains of step (3) ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value, realize the single absorbing boundary condition under 2D space coordinates.

2) as shown in Figure 2, for further improving the effect that suppresses edge reflection, the single absorption edge interlayer I in Fig. 1 in 2D absorbing boundary region 2 ₁j ₁k ₁l ₁arranged outside n-1 layer AWWE absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density under 2D space coordinates of the variable density under 2D space coordinates, by the mode of multi-ply linear weighting, calculate clathrum ω ₂=1,2 ..., the acoustic wavefield value on n, it specifically comprises the following steps:

(1) utilize the ACOUSTIC WAVE EQUATION computing grid layer ω of the variable density under 2D space coordinates ₂=0,1,2 ..., the acoustic wavefield value on n-1, and by the outermost layer absorption edge interlayer I in 2D absorbing boundary region 2 ₂j ₂k ₂l ₂(be clathrum ω ₂=n) the acoustic wavefield value on is preset as zero, by clathrum ω ₂=1,2 ..., n-1, the acoustic wavefield value on n is designated as

(2) according to the acoustic wavefield value on the outer boundary ABCD of 2D interior zone 1, utilize AWWE and the formula (7) of the variable density under 2D space coordinates, extrapolation calculates clathrum ω ₂acoustic wavefield value v on=1 ₁; Respectively according to the clathrum ω being obtained by step 1) ₂=1,2 ..., the acoustic wavefield value on n-1

utilize AWWE and the formula (7) of the variable density under 2D space coordinates, extrapolation calculates clathrum ω ₂=2 ..., n-1, the acoustic wavefield value on n

(3) the clathrum ω that will be obtained by step (1) ₂=1,2 ..., n-1, the acoustic wavefield value on n

with the clathrum ω being obtained by step (2) ₂=1,2 ..., n-1, the acoustic wavefield value on n

carry out linear weighted function, finally obtain clathrum ω after weighting ₂=1,2 ..., n-1, the acoustic wavefield value on n for:

$p_{{\mathrm{\&omega;}}_2}=\frac{{\mathrm{\&omega;}}_2}{n}{v}_{{\mathrm{\&omega;}}_2}+\frac{n-{\mathrm{\&omega;}}_2}{n}{u}_{{\mathrm{\&omega;}}_2},({\mathrm{\&omega;}}_2=\mathrm{1,2},...,n)---\left(12\right)$

3) as shown in Figure 3, based on step 1), 2D computer memory is expanded to 3D(three-dimensional) computer memory, in 3D computer memory, single absorption edge interlayer need to be set at the up, down, left, right, before and after six direction in 3D internal calculation region 11 simultaneously, to suppress the edge reflection from the sound wave of all directions; According to the acoustic wavefield value on the outer boundary ABCD of 3D interior zone 11, utilize 15 ° of one way wave equations amended under the AWWE treated side, 2D space coordinates of the variable density under 3d space coordinate system to process the array mode of 5 ° of one way wave equation processing angles under rib and 2D space coordinates, calculate the acoustic wavefield value on the single absorption edge interlayer of 3D interior zone 11 six directions.

Because the implementation method of the single absorption edge interlayer on the up, down, left, right, before and after six direction in 3D internal calculation region 11 is identical, therefore in 3D computer memory only to calculate the single absorption edge interlayer I of 3D interior zone 11 down directions ₃j ₃k ₃l ₃on acoustic wavefield value be that example describes, it specifically comprises the following steps:

(1), according to the acoustic wavefield value on the outer boundary ABCD in 3D internal calculation region 11, utilize the AWWE of the descending variable density under 3d space coordinate system to calculate additional clathrum

upper region M ₃n ₃o ₃p ₃four interior rib M ₃n ₃, interior rib N ₃o ₃, interior rib O ₃p ₃with interior rib P ₃m ₃on acoustic wavefield value.

The equation form of the AWWE of the descending variable density under 3d space coordinate system is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}+{\mathrm{<figure-callout\; id="1"\; label="cH"\; filenames="HDA0000451812910000042.png"\; state="\{\{state\}\}">cH</figure-callout>}}_1\frac{\&PartialD;}{\&PartialD;z})p={\mathrm{\&rho;c}}^2{H}_2(\frac{{\&PartialD;v}_x}{\&PartialD;x}+\frac{{\&PartialD;v}_y}{\&PartialD;y})\\ \frac{{\&PartialD;v}_x}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;x}\\ \frac{{\&PartialD;v}_y}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;y}\end{array}\right.---\left(13\right)$

In formula (13), x, y, z represents the coordinate under 3d space coordinate system, v _yrepresent the speed component of vibration in the y-direction that comprises auxiliary variable.

The finite difference computation scheme of the AWWE of the descending variable density under 3d space coordinate system is:

$\left\{\begin{array}{c}{p}_{i,j,k+\frac{1}{2}}^{n+1}+{\mathrm{\&alpha;}}_z{p}_{i,j,k+\frac{1}{2}}^{n+1}a=p_{i,j,k+\frac{1}{2}}^n+{\mathrm{\&alpha;}}_{z}a(p_{i,j,k}^{n+1}+p_{i,j,k}^n-p_{i,j,k+\frac{1}{2}}^n)\\ +\mathrm{\&rho;}{\mathrm{cH}}_2[{\mathrm{\&alpha;}}_x({\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_x\right)}_{i-\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}})+{\mathrm{\&alpha;}}_y({\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_y\right)}_{i,j-\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}})]\\ {\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},j,k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i+1,j,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\\ {\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i,j+\frac{1}{2},k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i,j+1,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\end{array}\right.---\left(14\right)$

In formula (14), x, y and z all represent the coordinate under 3d space coordinate system; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing;

v _xrepresent the speed component of vibration in the x-direction,

for calculating the auxiliary variable of introducing, v _yrepresent the speed component of vibration in the y-direction,

for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; H ₂for matrix, i, j and k are respectively the discrete grid block node under 3d space coordinate system, the discrete grid block node that n is the time; Matrix a is: a=(L ₁₁, L ₂₁, L ₃₁..., L _m1) ^t; α _x=c Δ t/ Δ x, α _y=c Δ t/ Δ y, α _z=c Δ t/ Δ z, Δ x, Δ y and Δ z are respectively the mesh spacing of x, y, z direction under 3d space coordinate system, the mesh spacing that Δ t is the time.

(2) at additional clathrum upper, according to four interior rib O ₃p ₃, M ₃n ₃, P ₃m ₃and N ₃o ₃on acoustic wavefield value, utilize respectively amendedly under 2D space coordinates move ahead, 15 ° of one way wave equations of rear row, left lateral and right lateral, calculate four outer rib V ₃w ₃, outer rib R ₃s ₃, outer rib X ₃z ₃with outer rib T ₃u ₃on acoustic wavefield value.

Under 2D space coordinates, amended moving ahead with the form of 15 ° of one way wave equations of rear row is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;y\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;\mathrm{<figure-callout\; id="2"\; label="t"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">t</figure-callout>}}^2}+\frac{c}{2}(\frac{1}{{\mathrm{<figure-callout\; id="2"\; label="c"\; filenames="HDA0000451812910000032.png"\; state="\{\{state\}\}">c</figure-callout>}}^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;y}^2})=0,---\left(15\right)$

In formula (15), amended moving ahead in 15 ° of one way wave equations

get "+", in 15 ° of one way wave equations of amended rear row

get "-".

(3) for suppressing angle point reflection, utilize 5 ° of one way wave equations under 2D space coordinates, respectively according to four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value on two adjoint points separately, calculates additional clathrum

upper four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value.

Calculate additional clathrum upper left front, left back, right back, right front four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value time 5 ° of one way wave equation forms under the 2D space coordinates utilized be respectively:

$\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;x}\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;y}+\frac{\sqrt{2}}{c}\frac{\&PartialD;p}{\&PartialD;t}=0---\left(16\right)$

(4) according to formula (7), respectively by the acoustic wavefield value substitution p on each rib in the outer boundary ABCD of 3D interior zone 11 _k, will add clathrum respectively

in acoustic wavefield value substitution on each rib the p that correspondence calculates _k+1value be the single absorption edge interlayer I of 3D interior zone 11 down directions ₃j ₃k ₃l ₃in acoustic wavefield value on each rib.

Utilize and the single absorption edge interlayer I that calculates 3D interior zone 11 down directions ₃j ₃k ₃l ₃on the identical processing mode of acoustic wavefield value, calculate the acoustic wavefield value on the single absorption edge interlayer of 3D interior zone 11 other five directions, realize the loading of the single absorption edge interlayer of AWWE of 3D variable density.

4) based on step 3), integrating step 2) implementation method on the AWWE blended absorbent border of 2D variable density in, the outside of the single absorption edge interlayer simultaneously arranging on the up, down, left, right, before and after six direction of 3D interior zone 11 all arranges multilayer absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density of the variable density under 3d space coordinate system, the mode by multi-ply linear weighting calculates clathrum ω ₃=1,2 ..., the acoustic wavefield value on n, the AWWE blended absorbent border of realizing 3D variable density.

Embodiment: Fig. 4 is that the combined strategy of individual layer AWWE and amended 15 ° of one way wave equations is for the assimilation effect schematic diagram of 2D uniform dielectric acoustic wavefield numerical simulation, Fig. 4 (a) and Fig. 4 (b) are phase acoustic wavefield sections in the same time in uniform dielectric, and show within the scope of identical colour code, all utilize single absorption edge interlayer.Fig. 4 (a) utilizes AWWE to process rib, gets 5 reference velocities and is equal to background velocity, utilizes 15 ° of one way wave equation processing angles simultaneously, and its implementation is shown in step 1); Fig. 4 (b) utilizes amended 15 ° of one way wave equations to process rib, 5 ° of one way wave equation processing angles.Clearly, the effect of the combination of AWWE and amended 15 ° of one way wave equations compacting edge reflection and angle point reflection is better than the combination of the middle low order one way wave equation of Fig. 4 (b), and in Fig. 4, white arrow indication is edge reflection.

Fig. 5 is that the AWWE blended absorbent border of 2D variable density is for the assimilation effect schematic diagram of uniform dielectric acoustic wavefield numerical simulation, wave field section in Fig. 5 (a) shows that wave field does not also arrive border, Fig. 5 (b) and 5(c) be mutually acoustic wavefield section in the same time, colour code shows in same range.Fig. 5 (b) utilizes individual layer AWWE absorbing boundary, and Fig. 5 (c) utilizes 10 layers of AWWE to form blended absorbent border, and its implementation is shown in step 2).The assimilation effect on the blended absorbent border of 10 layers of AWWE composition is obviously better than individual layer AWWE absorbing boundary.

Fig. 6 is Marmousi rate pattern, and corresponding inhomogeneous density model is generated by an experimental formula.

Fig. 7 is that the AWWE blended absorbent border of 2D variable density is for the assimilation effect schematic diagram of Marmousi model acoustic wavefield numerical simulation, Fig. 7 (a) and Fig. 7 (c) are two not acoustic wavefield in the same time sections, do not utilize absorbing boundary, therefore very obvious from the reflection of left margin and coboundary, in figure, white arrow indication is edge reflection; Fig. 7 (b) and Fig. 7 (d) correspond respectively to Fig. 7 (a) and Fig. 7 (c), utilize the AWWE of 10 layers of variable density to form blended absorbent border, and the effect of its compacting edge reflection is very remarkable.

Fig. 8 is the seismologic record receiving in the close a certain degree of depth on earth's surface in Fig. 6, Fig. 8 (a) and 8(b) corresponding to absorbing boundary and the situation of utilizing the AWWE blended absorbent border of 10 layers of variable density respectively.Contrast two figure, reconfirm the assimilation effect on AWWE blended absorbent border.

Fig. 9 be the AWWE blended absorbent border of 3D variable density for the assimilation effect schematic diagram of uniform dielectric acoustic wavefield numerical simulation, focus is positioned at upper surface center, the degree of depth is 20m.Fig. 9 (a) is the acoustic wavefield section of 155ms, shows that wave field does not also arrive left and right, front and back and lower boundary.Fig. 9 (b) and 9(c) be the acoustic wavefield section of 250ms, wherein Fig. 9 (b) utilizes individual layer AWWE absorbing boundary treated side, amended 15 ° of one way wave equations are processed rib, 5 ° of one way wave equation processing angles, obtained certain assimilation effect, but still visible faint reflection, as shown in white arrow in Fig. 9 (b); Fig. 9 (c) utilizes 10 layers of AWWE absorption edge interlayer to form blended absorbent border, and its assimilation effect is obviously better than individual layer AWWE absorbing boundary.Fig. 9 (d) and 9(e) be the acoustic wavefield section of 325ms, utilizes respectively individual layer AWWE and 10 layers of AWWE blended absorbent border, and the assimilation effect on 10 layers of AWWE blended absorbent border is better.

Figure 10 is 3D rate pattern heterogeneous, and density model heterogeneous is generated by an experimental formula.

Figure 11 is that the AWWE blended absorbent border of 3D variable density is for the assimilation effect schematic diagram of nonhomogeneous media acoustic wavefield numerical simulation, Figure 11 (a) and 11(c) difference corresponding diagram 11(b) and 11(d), Figure 11 (a) and 11(c) in without absorbing boundary, edge reflection is very obvious, as shown in the white arrow in figure (c); Figure 11 (b) and 11(d) utilize the AWWE blended absorbent border of 10 layers of variable density, its assimilation effect is very good.

The various embodiments described above are only for illustrating the present invention; wherein the structure of each parts, connected mode and method step etc. all can change to some extent; every equivalents of carrying out on the basis of technical solution of the present invention and improvement, all should not get rid of outside protection scope of the present invention.

Claims (10)

1. for an implementation method for the mixed absorbing boundary of variable density ACOUSTIC WAVE EQUATION, it comprises the following steps:

1) when the acoustic wavefield of simulation variable density, whole 2D computer memory is divided into 2D interior zone and comprises single absorption edge interlayer I ₁j ₁k ₁l ₁2D absorbing boundary region; Utilize the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on the outer boundary ABCD of 2D interior zone, according to the acoustic wavefield value on the outer boundary ABCD calculating, utilize the AWWE of the variable density under 2D space coordinates to process the array mode of 15 ° of one way wave equation processing angles amended under rib and 2D space coordinates, the single absorption edge interlayer I in calculating absorbing boundary region ₁j ₁k ₁l ₁on acoustic wavefield value;

2) the single absorption edge interlayer I in 2D absorbing boundary region ₁j ₁k ₁l ₁arranged outside n-1 layer AWWE absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density of the variable density under 2D space coordinates, the mode by multi-ply linear weighting calculates clathrum ω ₂=1,2 ..., the acoustic wavefield value on n;

3) 2D computer memory is expanded to 3D computer memory, on the up, down, left, right, before and after six direction in 3D internal calculation region, single absorbing boundary is set simultaneously, according to the acoustic wavefield value on the outer boundary ABCD of 3D interior zone, utilize amended 15 ° of one way wave equations under the AWWE treated side, 2D space coordinates of the variable density under 3d space coordinate system to process the array mode of 5 ° of one way wave equation processing angles under rib and 2D space coordinates, calculate the acoustic wavefield value on the single absorption edge interlayer of 3D interior zone six direction;

4) outside of the single absorption edge interlayer simultaneously arranging on the up, down, left, right, before and after six direction of 3D interior zone all arranges multilayer absorption edge interlayer, utilize the AWWE of ACOUSTIC WAVE EQUATION and the variable density of the variable density under 3d space coordinate system, the mode by multi-ply linear weighting calculates clathrum ω ₃=1,2 ..., the acoustic wavefield value on n, the AWWE blended absorbent border of realizing 3D variable density.

2. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 1, it is characterized in that: in described step 1), utilize the AWWE and amended 15 ° of modes that one way wave equation is combined of the variable density under 2D space coordinates, calculate the single absorption edge interlayer I in absorbing boundary region ₁j ₁k ₁l ₁on acoustic wavefield value, it comprises the following steps:

(1) the outer boundary ABCD of 2D interior zone comprises rib BC, rib CD, rib DA and rib AB, utilizes the ACOUSTIC WAVE EQUATION of the variable density under 2D space coordinates to calculate the acoustic wavefield value on each rib in outer boundary ABCD;

(2) according to the acoustic wavefield value on the rib BC of the outer boundary ABCD of 2D interior zone, rib CD, rib DA and rib AB, utilize respectively the AWWE of descending, right lateral under 2D space coordinates, up and left lateral variable density, calculate lower boundary J ₁k ₁, right margin K ₁l ₁, coboundary L ₁i ₁with left margin I ₁j ₁in four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value;

(3) utilize amended 15 ° of one way wave equations under 2D space coordinates, respectively according to additional clathrum E ₁f ₁g ₁h ₁the acoustic wavefield value of mid point M, N, O, P, Q, R, S, T, calculates additional clathrum E ₁f ₁g ₁h ₁in the outer end points W of four rib MN, OP, QR, ST ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value;

(4) the single absorption edge interlayer I obtaining according to step (2) ₁j ₁k ₁l ₁in four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value and the end points W that obtains of step (3) ₁, W ₂, W ₃, W ₄, W ₅, W ₆, W ₇, W ₈acoustic wavefield value, realize the AWWE absorbing boundary condition of the simple layer variable density under 2D space coordinates.

3. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 2, is characterized in that: in described step (2), and four rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on the computing method of acoustic wavefield value specifically comprise the following steps:

1. according to the acoustic wavefield value on rib BC, rib CD, rib DA and rib AB in the outer boundary ABCD of 2D interior zone, utilize the AWWE of descending, right lateral under 2D space coordinates, up and left lateral variable density, by finite difference, calculate additional clathrum E ₁f ₁g ₁h ₁acoustic wavefield value on middle rib MN, rib OP, rib QR, rib ST;

2. according to following formula

$p_{k+1}=2p_{k+\frac{1}{2}}-p_k,$

Respectively by the acoustic wavefield value substitution p on rib BC, rib CD, rib DA and rib AB _k, respectively by the acoustic wavefield value on rib MN, rib OP, rib QR and rib ST

the p that correspondence calculates _k+1value be single absorption edge interlayer I ₁j ₁k ₁l ₁middle rib M ₁n ₁, rib O ₁p ₁, rib Q ₁r ₁with rib S ₁t ₁on acoustic wavefield value.

4. a kind of implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 2 or claim 3, is characterized in that: the equation form of the AWWE of the descending and up variable density under described 2D space coordinates is:

$\left\{\begin{array}{c}(\frac{\&PartialD;}{\&PartialD;t}\&PlusMinus;{\mathrm{cH}}_1\frac{\&PartialD;}{\&PartialD;z})p={\mathrm{\&rho;c}}^2{H}_2\frac{\&PartialD;v_x}{\&PartialD;x}\\ \frac{\&PartialD;v_x}{\&PartialD;t}=\frac{1}{\mathrm{\&rho;}}\frac{\&PartialD;p}{\&PartialD;x}\end{array}\right.,$

In formula, in the equation of the AWWE of descending variable density

get "+", in the equation of the AWWE of up variable density

get "-"; X and z represent the coordinate under 2D space coordinates, and t represents time coordinate; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing; v _xrepresent the speed component of vibration in the x-direction, for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; H ₁and H ₂equal representing matrix;

X in the equation of the AWWE of descending variable density and z are exchanged, obtain the equation of the AWWE of right lateral variable density; X in the equation of the AWWE of up variable density and z are exchanged, obtain the equation of the AWWE of left lateral variable density.

5. a kind of implementation method of the mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 2 or claim 3, is characterized in that: under described 2D space coordinates, the form of amended left lateral and 15 ° of one way wave equations of right lateral is:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;x\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}+\frac{c}{2}(\frac{1}{c^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;x}^2})=0,$

In formula, in 15 ° of one way wave equations of amended left lateral

get "+", in 15 ° of one way wave equations of amended right lateral

get "-";

Change the x in 15 ° of one way wave equations of amended left lateral into z, obtain amended up 15 ° of one way wave equations under 2D space coordinates; Change the x in 15 ° of one way wave equations of amended right lateral into z, obtain amended descending 15 ° of one way wave equations under 2D space coordinates.

6. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as described in claim 1 or 2 or 3, is characterized in that: described step 2), and computing grid layer ω ₂=1,2 ..., the acoustic wavefield value on n, it specifically comprises the following steps:

(1) utilize the ACOUSTIC WAVE EQUATION computing grid layer ω of the variable density under 2D space coordinates ₂=1,2 ..., the acoustic wavefield value on n-1, and by clathrum ω ₂acoustic wavefield value on=n is preset as zero, clathrum ω ₂=1,2 ..., n-1, the acoustic wavefield value on n is designated as

(2) according to the acoustic wavefield value on the outer boundary ABCD of 2D interior zone, utilize the variable density under 2D space coordinates AWWE and

extrapolation calculates clathrum ω ₂acoustic wavefield value v on=1 ₁; Respectively according to the clathrum ω being obtained by step 1) ₂=1,2 ..., the acoustic wavefield value on n-1 utilize the variable density under 2D space coordinates AWWE and

extrapolation calculates clathrum ω ₂=2 ..., n-1, the acoustic wavefield value on n

(3) the clathrum ω that will be obtained by step (1) ₂=1,2 ..., n-1, the acoustic wavefield value on n

with the clathrum ω being obtained by step (2) ₂=1,2 ..., n-1, the acoustic wavefield value on n

carry out linear weighted function, obtain clathrum ω after weighting ₂=1,2 ..., n-1, the acoustic wavefield value on n

for:

$p_{{\mathrm{\&omega;}}_2}=\frac{{\mathrm{\&omega;}}_2}{n}{v}_{{\mathrm{\&omega;}}_2}+\frac{n-{\mathrm{\&omega;}}_2}{n}{u}_{{\mathrm{\&omega;}}_2},({\mathrm{\&omega;}}_2=\mathrm{1,2},...,n).$

7. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as described in claim 1 or 2 or 3, it is characterized in that: in described step 3), the implementation method of the single absorption edge interlayer on the up, down, left, right, before and after six direction in 3D internal calculation region is identical, wherein calculates the single absorption edge interlayer I of 3D interior zone down direction ₃j ₃k ₃l ₃on acoustic wavefield value comprise the following steps:

(1) according to the acoustic wavefield value on the outer boundary ABCD in 3D internal calculation region, utilize the AWWE of the descending variable density under 3d space coordinate system, calculate additional clathrum

upper region M ₃n ₃o ₃p ₃four interior rib M ₃n ₃, interior rib N ₃o ₃, interior rib O ₃p ₃with interior rib P ₃m ₃on acoustic wavefield value;

(2) according to four interior rib M ₃n ₃, interior rib N ₃o ₃, interior rib O ₃p ₃with interior rib P ₃m ₃on acoustic wavefield value, utilize respectively amendedly under 2D space coordinates move ahead, 15 ° of one way wave equations of rear row, left lateral and right lateral, calculate additional clathrum

upper four outer rib V ₃w ₃, outer rib R ₃s ₃, outer rib X ₃z ₃with outer rib T ₃u ₃on acoustic wavefield value;

(3) utilize 5 ° of one way wave equations under 2D space coordinates, respectively according to four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value on two adjoint points separately, calculates additional clathrum

upper four angle point E ₃, F ₃, G ₃, H ₃acoustic wavefield value;

(4) according to formula

respectively by the acoustic wavefield value substitution p on each rib in the outer boundary ABCD of 3D interior zone 11 _k, will add clathrum respectively in acoustic wavefield value substitution on each rib

the p that correspondence calculates _k+1value be the single absorption edge interlayer I of 3D interior zone down direction ₃j ₃k ₃l ₃in acoustic wavefield value on each rib.

8. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 7, is characterized in that: in described step (1), the finite difference computation scheme of the AWWE of the descending variable density under 3d space coordinate system is:

$\left\{\begin{array}{c}{p}_{i,j,k+\frac{1}{2}}^{n+1}+{\mathrm{\&alpha;}}_z{p}_{i,j,k+\frac{1}{2}}^{n+1}a=p_{i,j,k+\frac{1}{2}}^n+{\mathrm{\&alpha;}}_{z}a(p_{i,j,k}^{n+1}+p_{i,j,k}^n-p_{i,j,k+\frac{1}{2}}^n)\\ +\mathrm{\&rho;}{\mathrm{cH}}_2[{\mathrm{\&alpha;}}_x({\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_x\right)}_{i-\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}})+{\mathrm{\&alpha;}}_y({\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}-{\left(v_y\right)}_{i,j-\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}})]\\ {\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_x\right)}_{i+\frac{1}{2},j,k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},j,k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i+1,j,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\\ {\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n+\frac{1}{2}}={\left(v_y\right)}_{i,j+\frac{1}{2},k+\frac{1}{2}}^{n-\frac{1}{2}}+\frac{1}{{\mathrm{\&rho;}}_{i+\frac{1}{2},j,k+\frac{1}{2}}}\frac{\mathrm{\&Delta;t}}{\mathrm{\&Delta;x}}(p_{i,j+1,k+\frac{1}{2}}^n-p_{i,j,k+\frac{1}{2}}^n)\end{array}\right.,$

In formula, x, y and z all represent the coordinate under 3d space coordinate system; P represents the m dimensional vector that comprises auxiliary variable, p=(p, p ₁... p _m-1) ^t, p represents sonic pressure field, (p ₁, p ₂..., p _m-1) for calculating the auxiliary variable of introducing;

v _xrepresent the speed component of vibration in the x-direction,

for calculating the auxiliary variable of introducing, v _yrepresent the speed component of vibration in the y-direction,

for calculating the auxiliary variable of introducing; ρ represents density, and c represents background velocity; H ₂for matrix, i, j and k are respectively the discrete grid block node under 3d space coordinate system, the discrete grid block node that n is the time; Matrix a is: a=(L ₁₁, L ₂₁, L ₃₁..., L _m1) ^t; α _x=c Δ t/ Δ x, α _y=c Δ t/ Δ y, α _z=c Δ t/ Δ z, Δ x, Δ y and Δ z be under 3d space coordinate system mesh spacing, the mesh spacing that Δ t is the time.

9. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 7, is characterized in that: in described step (2), and amended move ahead and the form of 15 ° of one way wave equations of rear row is under 2D space coordinates:

$\&PlusMinus;\frac{{\&PartialD;}^{2}p}{\&PartialD;y\&PartialD;t}-\frac{1}{c}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}+\frac{c}{2}(\frac{1}{c^2}\frac{{\&PartialD;}^{2}p}{{\&PartialD;t}^2}-\frac{{\&PartialD;}^{2}p}{{\&PartialD;y}^2})=0;$

In formula, amended moving ahead in 15 ° of one way wave equations

get "+", in 15 ° of one way wave equations of amended rear row get "-".

10. the implementation method of a kind of mixed absorbing boundary for variable density ACOUSTIC WAVE EQUATION as claimed in claim 7, is characterized in that: in described step (3), 5 ° of one way wave equation forms under 2D space coordinates are:

$\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;x}\&PlusMinus;\frac{\&PartialD;p}{\&PartialD;y}+\frac{\sqrt{2}}{c}\frac{\&PartialD;p}{\&PartialD;t}=0.$

Images