Seismic wave forward modeling and two-way wave illumination method of two-phase viscoelastic TTI (Transmission Time Interval) medium

Abstract

The invention discloses a seismic wave forward modeling simulation and two-way wave illumination method for a two-phase viscoelastic TTI medium, and relates to the technical field of geophysical exploration of petroleum, and the method comprises the steps: introducing a Biot two-phase medium theory to describe a two-phase medium; deriving a constitutive equation of the TTI medium; obtaining a stress-strain relationship of the solid skeleton; deriving a wave equation suitable for the biphase viscoelastic TTI medium; carrying out seismic wave forward modeling simulation on the derived first-order velocity stress equation of the biphase viscoelastic TTI medium by adopting a high-order staggered grid finite difference method; establishing an illumination operator, defining a source illumination operator and a reverse illumination operator, and calculating source illuminance through forward modeling; based on a reciprocal principle, a reverse illumination result is analyzed, and an excitation position and shot point density are optimized. According to the method, the TTI medium and the Biot biphase viscoelastic medium are combined for the first time, effective exploration illumination of complex biphase media with viscoelasticity and TTI anisotropy such as oil-gas-containing stratums is achieved, and the defect that currently, propagation characteristics of seismic waves in the media are insufficiently recognized is overcome.

Description

Technical Field

[0001] The invention relates to the technical field of petroleum geophysical exploration, and in particular to a seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium. Background Art

[0002] In oil and gas exploration, oil and gas-bearing strata, due to their unique stratification and directional development of fractures, coexist with solid and fluid phases, forming a two-phase complex medium with viscoelasticity and inclined transverse anisotropy (TTI).

[0003] Therefore, for the needs of oil and gas exploration, it is particularly important to understand the propagation characteristics of seismic waves in this complex medium, and it is urgent to study a seismic wave forward simulation and two-way wave illumination method for two-phase viscoelastic TTI medium. Summary of the invention

[0004] In order to achieve forward simulation of complex oil and gas media and describe underground media more realistically and accurately, the present invention discloses a seismic wave forward simulation and two-way wave illumination method for a two-phase viscoelastic TTI medium.

[0005] To achieve the above object, the present invention adopts the following technical solutions:

[0006] A seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium comprises the following steps:

[0007] Step S1, introducing the Biot two-phase medium theory to describe the two-phase medium, and using the coupling coefficient to represent the interaction between the volume change of the solid phase and the liquid phase;

[0008] Step S2, considering the anisotropy of the two-phase medium, taking the TTI medium as the solid skeleton, and deriving the constitutive equation of the TTI medium;

[0009] Step S3, introduce Carcione's viscoelastic anisotropic medium model and combine it with the Boltzmann superposition principle to obtain the stress-strain relationship of the solid skeleton:

[0010] Step S4, embedding the viscoelastic constitutive relation of the solid phase into Biot's two-phase medium theory, and deriving a wave equation applicable to the two-phase viscoelastic TTI medium;

[0011] Step S5, using a high-order staggered grid finite difference method to perform seismic wave forward simulation on the derived first-order velocity stress equation of the two-phase viscoelastic TTI medium;

[0012] Step S6, establishing an illumination operator, defining source direction and reverse illumination operators, calculating source direction illumination through forward simulation, and analyzing the illumination effects on the target area under different configurations;

[0013] Step S7: Based on the reciprocity principle, analyze the reverse illumination results, identify the surface locations that contribute more to the target area, and optimize the excitation positions and shot point density to achieve optimal coverage.

[0014] Optionally, in step S2, the constitutive equation of the TTI medium is:

[0015] ;

[0016] In the formula, is the TTI medium stress vector, is the TTI medium strain vector.

[0017] Optionally, in step S3, the stress-strain relationship of the solid skeleton is:

[0018] ;

[0019] In the formula, is the stress vector, is the first-order time derivative of the strain vector, represents the relaxation matrix;

[0020] Separate the integral and non-integral terms in the equation and substitute them into the equation ,get:

[0021] ;

[0022] In the formula, represents the elasticity matrix of the non-integral term, represents the viscoelastic matrix of the integral term, Represents a memory variable, They represent the expansion and contraction motion and shear motion of particles in the medium, respectively. Represents the number of units that make up the generalized standard linear device.

[0023] Optionally, in step S4, the wave equation applicable to the two-phase viscoelastic TTI medium is:

[0024] ;

[0025] In the formula, is the displacement vector; express The second derivative with respect to time; is the body force vector, where , , represents the solid phase fraction, , , represents the liquid phase fraction, Compared to A 7×7 matrix with a row vector and a column vector with 0 elements is added respectively; is the density matrix; is a partial differential operator; is the dissipation coefficient matrix; is the sum and difference operator.

[0026] Optionally, in step S5, the first-order time derivatives of the stress and strain of the wave equation of the two-phase viscoelastic TTI medium are obtained to obtain the first-order velocity stress equation of the two-phase viscoelastic TTI medium; in the two-dimensional case, the velocity stress equation is:

[0027] ;

[0028] ;

[0029] In the formula, , , ; is the solid effective mass of the medium volume element; is the mass coupling coefficient between the solid phase and the fluid phase, also called coupling mass; is the effective mass of the liquid phase per volume element of the medium; and Solid phase Direction and Direction speed; represents the dissipation coefficient; and Liquid phase Direction and Direction speed; ;

[0030] The first-order time derivative equation of the memory variable is:

[0031] ;

[0032] In the formula, relaxation time It can be calculated by the quality factor, the calculation formula is:

[0033] ;

[0034] In the formula, , and are the longitudinal wave quality factor and the shear wave quality factor, It is the reciprocal of the main frequency of the source wavelet.

[0035] Optionally, in step S6, the energy and illumination are calculated based on the two-way wave forward modeling, and the wave field value is recorded as , mark the shadow area of ​​the observation system source illumination (also called forward illumination) result as the pseudo reverse excitation point , the area where the surface observation system is located is recorded as , The internal position coordinates are generally expressed as , the source excitation point of source-directed illumination is recorded as , then the illumination of the observation system to the reverse excitation point is:

[0036] ;

[0037] The illumination function of the reversely excited wave field on the surface observation system is:

[0038] ;

[0039] According to the coordinate point corresponding to the function extreme value of the above formula , by the function Instructions, coordinate points The backlighting energy at is the maximum in the surface observation system.

[0040] Optionally, in step S7, the encrypted surface observation system source excitation point is recorded as , underground medium model area The coordinates of any point inside are marked as , then the illumination function of the optimized observation system for the medium model is:

[0041] ;

[0042] The above formula is the optimized lighting function for lighting enhancement in shadow areas after two-way wave lighting analysis, i.e., the lighting operator.

[0043] According to the derivation logic of the equilibrium equation of motion of two-phase medium, it is not difficult to see that the energy of two-phase medium is composed of the energy of the solid phase and liquid phase, as well as the energy lost by the solid-liquid coupling effect. If the solid-phase wave field of two-phase medium is recorded as Quantity and Quantity and , liquid phase wave field Quantity and Quantity and , then the solid phase energy is:

[0044] ;

[0045] The liquid phase energy is:

[0046] ;

[0047] The energy of the coupling effect is:

[0048] ;

[0049] So the total energy of the two-phase medium is:

[0050] ;

[0051] In two-dimensional illumination, if the coordinates of each position in the surface observation range are generally expressed as , the excitation positions of the observation system are , in the two-phase viscoelastic medium model region The coordinates of any point are , then according to the principle of statistical illumination method, the source-directed illumination function of the observation system for the two-phase medium is:

[0052] ;

[0053] It can be used as the source-direction illumination operator of the two-phase viscoelastic TTI medium of the present invention;

[0054] The target point of the source illumination in the shadow area is , the shaded area is , several excitation points are arranged in the shadow area for reverse lighting, and the lighting function within the observation range of the surface obtained by statistics is:

[0055] ;

[0056] It can be used as a reverse illumination operator for the two-phase viscoelastic TTI medium of the present invention;

[0057] remember The coordinates of the extreme points are , according to the lighting optimization scheme based on the reciprocity principle, The seismic waves excited at the location contribute the most to the illumination of the shadow area, so the original observation system's shot point excitation position is adjusted. The shot points are encrypted and the observation system is optimized. The encrypted surface observation system is recorded as The illumination function of the optimized and encrypted observation system for the two-phase medium is:

[0058] .

[0059] The beneficial effect of the present invention is that the present invention is based on the Biot two-phase medium theory, introduces the Carcione system and combines the constitutive relationship of the TTI medium to describe the viscoelastic anisotropy of the solid skeleton. By embedding the constitutive relationship of the solid phase into the Biot two-phase system, the wave equation of the two-phase viscoelastic TTI medium is derived; then, the first-order velocity-stress equation of the two-phase viscoelastic TTI medium is simulated by seismic wave forward modeling using a high-order staggered grid finite difference method, and a two-phase illumination operator is established using the derived wave equation, thereby performing two-way wave illumination on the two-phase viscoelastic TTI medium; finally, the wave field characteristics of the two-phase viscoelastic TTI medium are analyzed through numerical simulation results. The present invention realizes effective exploration illumination of complex two-phase media with viscoelasticity and TTI anisotropy, such as oil and gas formations, and makes up for the current lack of understanding of the propagation characteristics of seismic waves in such media. BRIEF DESCRIPTION OF THE DRAWINGS

[0060] Figure 1 It is a flow chart of a seismic wave forward modeling and two-way wave illumination method of a two-phase viscoelastic TTI medium according to the present invention;

[0061] Figure 2 It is the wave field degradation diagram of the two-phase viscoelastic TTI medium;

[0062] Figure 3 The wave field degradation diagram of the two-phase viscoelastic TTI medium after adjusting the order;

[0063] Figure 4 Establish the parameters of the two-phase viscoelastic TTI medium model for the depression model construction of a certain work area;

[0064] Figure 5A This is the result diagram of the solid phase operator source illumination of the depression model in this work area;

[0065] Figure 5B This is the result diagram of the source direction illumination of the liquid phase operator of the depression model in this work area;

[0066] Figure 5C This is the source-direction illumination result diagram of the two-phase operator of the depression model in this work area;

[0067] Fig. 6A This is the source-direction illumination diagram of the depression model in the work area;

[0068] Figure 6B This is the reverse lighting diagram of the depression model in the work area;

[0069] Figure 7 This is the normalized contour map of the energy below 400 m depth from the source direction of the depression model in the work area;

[0070] Figure 8 Back-illuminate the surface energy curve for the depression model of the work area;

[0071] Fig. 9 The lighting results after optimizing and encrypting the sunken model of the work area;

[0072] Fig.10 The energy normalization below 400 m depth before and after the sag model of the work area is densified;

[0073] Fig.11 The distribution of the parameter values ​​of the Marmours model of two-phase viscoelastic TTI properties;

[0074] Fig. 12A The source-directed illumination result of the Marmours model of the two-phase viscoelastic TTI property obtained by the method of the present invention;

[0075] Fig. 12B This is the energy distribution map of the Marmours model below 1200 m depth;

[0076] Fig.13A The reverse illumination result of the Marmours model of the two-phase viscoelastic TTI property obtained by the method of the present invention;

[0077] Fig. 13B is the normalized curve of reverse illumination surface energy;

[0078] Fig.14A This is the source illumination result of the Marmours model after shot point optimization and encryption;

[0079] Fig. 14B Energy normalization for the part below 1200 m depth after encryption of the Marmours model;

[0080] Fig.15 Comparison of normalized energy before and after optimization of the Marmous model below 1200 m depth. DETAILED DESCRIPTION

[0081] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, rather than all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention. Therefore, the following detailed description of the embodiments of the present invention provided in the drawings is not intended to limit the scope of the invention claimed for protection, but merely represents the selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0082] The present invention proposes a seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium. Figure 1 As shown, the specific steps include:

[0083] Step S1, introducing the Biot two-phase medium theory to describe the two-phase medium; in the Biot theory, the geometric equation of the two-phase medium (i.e., the Cauchy equation) is:

[0084] (1);

[0085] In the formula, The solid phase displacement of the medium is expressed by represents the liquid phase displacement of the medium, e represents the solid phase strain, ε represents the liquid phase strain, ;

[0086] The general form of the constitutive equation for a two-phase anisotropic medium is:

[0087] (2);

[0088] In the formula, represents the solid elastic matrix, and the matrix elements are =1,2,…6; and are the solid and liquid phase force vectors, respectively, and are the solid and liquid phase deformation vectors respectively; represents the coupling coefficient vector, represents the liquid phase pressure coefficient vector; =1,2,…,6 is the coupling coefficient, which represents the interaction between the volume changes of the solid and liquid phases; R is the liquid phase pressure coefficient, which represents a measure of the force required on the fluid to force a certain volume of fluid to flow into the volume element of the two-phase medium while keeping the volume constant.

[0089] Step S2: Taking the TTI medium as the solid skeleton, deriving the constitutive equation of the TTI medium.

[0090] Based on the anisotropy of the two-phase medium, taking the TTI medium as the solid skeleton, the constitutive equation of the TTI medium is first introduced. The derivation process of the constitutive equation of the TTI medium is as follows:

[0091] Assume that the constitutive equation of VTI medium is:

[0092] (3);

[0093] In the formula, is the VTI medium stress vector, is the VTI medium strain vector, is the elastic matrix of the VTI medium; according to the Bond transformation, the elastic matrix is ​​transformed from the constitutive coordinate system to the observation coordinate system to obtain the elastic matrix of the TTI medium:

[0094] (4);

[0095] In the formula, is the polarization angle of the TI medium; is the azimuth angle of the TI medium; is the polarization anisotropy matrix, is the azimuthal anisotropy matrix;

[0096] Substituting equation (4) into equation (3), the constitutive equation of TTI medium is obtained as follows:

[0097] (5);

[0098] In the formula, is the TTI medium stress vector, is the TTI medium strain vector.

[0099] Step S3: Carcione's viscoelastic anisotropic medium model is introduced, and combined with the Boltzmann superposition principle, the stress-strain relationship of the solid skeleton is obtained.

[0100] Carcione introduced the superposition principle of the cumulative integral method and the anisotropy of the medium to obtain a viscoelastic anisotropic medium model based on a generalized standard linear body. According to Carcione's viscoelastic anisotropy theory and Boltzmann's superposition principle, the stress-strain relationship of the solid skeleton is:

[0101] (6);

[0102] In the formula, is the stress vector, is the first-order time derivative of the strain vector, represents the relaxation matrix.

[0103] Step S4: embed the viscoelastic constitutive relation of the solid phase into Biot's two-phase medium theory to derive the wave equation applicable to the two-phase viscoelastic TTI medium, which specifically includes:

[0104] Separating the integral term and the non-integral term in equation (6) and substituting them into equation (5), we obtain:

[0105] (7);

[0106] In the formula, represents the elasticity matrix of the non-integral term, represents the viscoelastic matrix of the integral term, Represents a memory variable, They represent the expansion and contraction motion and shear motion of particles in the medium, respectively. Indicates the number of units that make up the generalized standard linear device;

[0107] By inserting the viscoelastic TTI solid skeleton into the mechanism of the two-phase medium, ignoring the viscosity of the fluid, and replacing the solid term in equation (2) with equation (7), the constitutive equation of the two-phase viscoelastic TTI medium based on the standard linear viscoelastic body is obtained:

[0108] (8);

[0109] In the formula, , =1,…, , without considering the fluid viscosity, we have ;

[0110] Using the principles of Newtonian mechanics to analyze the dynamic characteristics of a two-phase medium, take any volume element of the medium, and the kinetic energy of the unit cube is:

[0111] (9);

[0112] In the formula, is the solid effective mass of the medium volume element; is the mass coupling coefficient between the solid phase and the fluid phase, also called coupling mass; is the effective mass of the liquid phase per volume element of the medium;

[0113] Due to the interaction between the solid phase and the liquid phase, there is energy loss per unit volume of the two-phase medium relative to the single-phase medium, so the dissipation function is defined as:

[0114] (10);

[0115] Where B is the dissipation coefficient, satisfying ; is the fluid viscosity coefficient, is the permeability;

[0116] After considering the energy dissipation of the solid and liquid phases, the energy of the volume element is:

[0117] (11);

[0118] Momentum of the solid phase for:

[0119] (12);

[0120] Momentum of the liquid phase for:

[0121] (13);

[0122] The force on the solid phase in the x direction is:

[0123] (14);

[0124] The force on the liquid phase in the x direction is:

[0125] (15);

[0126] According to the relationship between kinetic energy and momentum, the relationship between momentum and force, and equation (11), the force in the solid phase in the x direction is obtained as:

[0127] (16);

[0128] The force on the liquid phase in the x direction is:

[0129] (17);

[0130] Considering the anisotropy of the medium, the motion of the solid and liquid phases satisfies the generalized Darcy's law, and the dissipation coefficient is a symmetric matrix B, expressed as:

[0131] = (18);

[0132] The stress-displacement relationship does not consider viscoelasticity. The simultaneous equations of stress and momentum are used to express the motion equilibrium equation of the two-phase viscoelastic TTI medium. The stress-displacement equation of the solid phase is:

[0133] (19);

[0134] The stress-displacement equation of the liquid phase is:

[0135] (20);

[0136] Substituting equation (1) and constitutive equation (8) into the kinematic equilibrium equations (19) and (20), the wave equation of the two-phase viscoelastic TTI medium expressed in displacement form is obtained:

[0137] (twenty one);

[0138] In the formula, is the displacement vector; express The second derivative with respect to time; is the body force vector, where , , represents the solid phase fraction, , , represents the liquid phase fraction, Compared to A 7×7 matrix with a row vector and a column vector with 0 elements is added respectively; is the density matrix; is a partial differential operator; is the dissipation coefficient matrix; is a sum and difference operator;

[0139] In formula (21), the correlation coefficient, operator matrix and calculation formula are:

[0140] ;

[0141] ;

[0142] ;

[0143] ;

[0144] ;

[0145] ;

[0146] ;

[0147] In the formula, , , , , , , , ; represents the number of viscoelastic units in the solid phase, represents the number of viscoelastic units in the liquid phase. The values ​​of M and N determine the relaxation function and viscoelastic properties in the solid and liquid phases, respectively. , =1,…,6 represents the elastic coefficient, represents the relaxation function, and They correspond to expansion and contraction strain and shear strain, that is, the relaxation functions of longitudinal waves and transverse waves, respectively. and They represent the standard linear solids The strain relaxation time and stress relaxation time of the device assembly; is the inverse of the center frequency, is the longitudinal or transverse wave quality factor corresponding to the center frequency, Corresponding to the number of device assembly units in the viscoelastic model.

[0148] Step S5, using a high-order staggered grid finite difference method to perform seismic wave forward simulation on the derived first-order velocity stress equation of the two-phase viscoelastic TTI medium; taking the first-order time derivatives of the stress and strain of the wave equation of the two-phase viscoelastic TTI medium to obtain the first-order velocity stress equation of the two-phase viscoelastic TTI medium.

[0149] In the two-dimensional case, the velocity stress equation is:

[0150] (twenty two);

[0151] (twenty three);

[0152] In formula (22), , , ; is the solid effective mass of the medium volume element; is the mass coupling coefficient between the solid phase and the fluid phase, also called coupling mass; is the effective mass of the liquid phase per volume element of the medium; and Solid phase Direction and Direction speed; represents the dissipation coefficient; and Liquid phase Direction and The speed in the direction; in formula (23), ;

[0153] The first-order time derivative equation of the memory variable is:

[0154] (twenty four);

[0155] In the formula, relaxation time Calculated by quality factor, the calculation formula is:

[0156] (25);

[0157] In the formula, , and are the longitudinal wave quality factor and the shear wave quality factor, is the reciprocal of the main frequency of the source wavelet;

[0158] The staggered grid finite difference method is used to differentiate equations (22), (23) and (24) to obtain iterative formulas for calculating velocity, stress and memory variables. The perfectly matched layer (PML) method is used in the forward modeling to decompose the velocity and stress into two components parallel and perpendicular to the artificial boundary, and substituted them into the difference equation. By introducing the attenuation factor, the seismic wave is strongly attenuated in the artificially established attenuation layer to eliminate boundary reflections.

[0159] In order to make the difference stable, relevant stability conditions need to be met. Taking into account the viscoelasticity and TI anisotropy of the medium, the necessary stability conditions for the forward simulation of the first-order velocity stress equation of the two-phase viscoelastic TTI medium using the high-order staggered grid finite difference method are:

[0160] (26);

[0161] In the formula, is the time sampling interval, , are the spatial steps in the x-direction and z-direction, i.e., the grid spacing, is the 2N-order spatial staggered grid finite difference coefficient;

[0162] Using the Routh criterion, the stability condition that needs to be satisfied by the staggered grid finite difference forward modeling of viscoelastic media is derived as follows:

[0163] (27);

[0164] By analyzing the stability of the staggered grid finite difference method for viscoelastic TTI media, the stability condition is obtained as follows:

[0165] (28);

[0166] In the formula, .

[0167] Step S6: Establish an illumination operator, define source direction and reverse illumination operators, calculate source direction illumination by forward simulation, and analyze the illumination effect on the target area under different configurations; use the introduced two-phase viscoelastic TTI medium wave equation to perform illumination work, the process is as follows:

[0168] Assume that seismic waves propagate in underground space, a certain area underground Any particle is denoted as , the wave field of the particle is recorded as , the observation system has excitation points, each source is represented by the symbol ,…, According to the relationship between seismic wave energy and wave field value, particle exist The energy at that moment is recorded as , assuming that the total number of time sampling points in the whole process of seismic wave propagation in the forward simulation is , then a single source The illumination of the particle is:

[0169] (29);

[0170] The illumination of the particle by the entire observation system is:

[0171] (30);

[0172] Statistical observation system for model area The illumination of all particles in the medium is obtained by the observation system:

[0173] (31);

[0174] Equation (31) is the expression of the illumination function, also known as the illumination operator of the observation system on the medium model. The illumination process of the observation system excited from the surface to the underground medium is called source-directed illumination.

[0175] Based on the two-way wave forward modeling, the energy and illumination are calculated, and the wave field value is recorded as , mark the shadow area of ​​the observation system source illumination (also called forward illumination) result as the pseudo reverse excitation point , the area where the surface observation system is located is recorded as , The internal position coordinates are expressed as , the source excitation point of source-directed illumination is recorded as , then the illumination of the observation system to the reverse excitation point is:

[0176] (32);

[0177] The illumination function of the reversely excited wave field on the surface observation system is:

[0178] (33);

[0179] According to the coordinate point corresponding to the function extreme value of formula (33) , by the function Indication, coordinate point The backlighting energy at is the maximum in the surface observation system.

[0180] Step S7: Based on the reciprocity principle, analyze the reverse illumination results, identify the surface positions that contribute more to the target area, and optimize the excitation positions and shot point density to achieve the best coverage. Based on the reciprocity principle followed by seismic wave propagation, reverse illumination is performed with the target area as the source. The energy results on the surface can indicate the most favorable excitation position for improving the wave field energy in the target area, thereby guiding the adjustment of the excitation position and the number of shot points of the observation system, thereby optimizing the observation system.

[0181] According to the inverse principle, the seismic wave excited by the coordinate point with the extreme value in equation (33) affects the target point The lighting energy is the largest, that is, if The density of shot points can enhance the illumination contribution to the shadow area. The densified surface observation system source excitation point is , underground medium model area The coordinates of any point inside are marked as , then the illumination function of the optimized observation system for the medium model is:

[0182] (34);

[0183] Formula (34) is the optimized lighting function for lighting enhancement in shadow areas after two-way wave lighting analysis, i.e., the lighting operator;

[0184] The energy of the two-phase medium is composed of the energy of the solid phase and the liquid phase, as well as the energy lost due to the solid-liquid coupling effect. The x component and z component of the solid phase wave field of the two-phase medium are denoted as and , the x-component and z-component of the liquid phase wave field are and , the solid phase energy is:

[0185] (35);

[0186] The liquid phase energy is:

[0187] (36);

[0188] The energy of the coupling effect is:

[0189] (37);

[0190] The total energy of the two-phase medium is:

[0191] (38);

[0192] In two-dimensional illumination, the coordinates of each position in the surface observation range are recorded as , the excitation positions of the observation system are , in the two-phase viscoelastic medium model region The coordinates of any point are , then according to the principle of statistical illumination method, the source-directed illumination function of the observation system for the two-phase medium is:

[0193] (39);

[0194] In formula (39), As a source-directed illumination operator of the two-phase viscoelastic TTI medium of the present invention;

[0195] The target point of the source illumination in the shadow area is , the shaded area is , several excitation points are arranged in the shadow area for reverse lighting, and the lighting function within the observation range of the surface obtained by statistics is:

[0196] (40);

[0197] In formula (40), As a reverse illumination operator of the two-phase viscoelastic TTI medium of the present invention;

[0198] remember The coordinates of the extreme points are , according to the lighting optimization scheme based on the reciprocity principle, The seismic waves excited at the location contribute the most to the illumination of the shadow area, so the original observation system's shot point excitation position is adjusted. The shot points are encrypted and the observation system is optimized. The encrypted surface observation system is recorded as ; The illumination function of the optimized and encrypted observation system for the two-phase medium is:

[0199] (41);

[0200] Adjust the excitation position and shot point density in the observation system, especially in the shadow area where the lighting effect is originally weak. Through the illumination analysis of reverse excitation, select the energy maximum point for encryption, optimize the design of seismic exploration, and ensure the best coverage of seismic wave energy in the target area.

[0201] Application Examples

[0202] (a) Wave field numerical simulation analysis based on parameter changes

[0203] Create 301 A homogeneous medium model with a grid size of 301, using a Ricker wavelet longitudinal wave source at the center of the model The source frequency is 25 Hz, the time sampling interval is 0.5 ms, and the space sampling interval is = =10 m. The verification work starts with the forward modeling of the two-phase viscoelastic medium model. After that, each subsequent stage of the forward modeling sets the relevant parameters to zero based on the previous step. The wave equation derived in this application example is verified and analyzed based on the forward modeling results of each step. The initial parameters of the two-phase viscoelastic TTI medium model are shown in Table 1.

[0204] Table 1 Parameters of the uniform model of two-phase viscoelastic TTI medium

[0205]

[0206] Note: , The unit is m / s; , , The unit is kg / m 3 ; , The unit is 10 9 kg·m -1 ·s -2 ; The unit is kg·m -3 ·s -1 ; The remaining parameters are dimensionless.

[0207] like Figure 2 As shown in the figure, the wave field snapshots of each step in the process of the two-phase viscoelastic TTI medium first degenerating into various types and simple two-phase media, and finally degenerating into a single-phase uniform isotropic fully elastic medium. The specific steps are: from left to right, first set the polarization angle θ=0, the energy of the wave field of the two-phase medium does not change much, the wavefront morphology of the longitudinal and transverse waves changes, the symmetry axis rotates to the direction parallel to the z-axis, and it has the anisotropy and transverse wave splitting characteristics of fast longitudinal waves, and there are also slow longitudinal waves. Its wave field characteristics are consistent with those of the two-phase viscoelastic VTI medium, indicating that when the polarization angle is zero, the equation degenerates into the wave equation of the two-phase viscoelastic VTI medium; secondly, on the basis of θ=0, set the anisotropy parameter = = =0, the total energy of the wave field still does not change much, but the wave front shape has changed greatly. The fast longitudinal wave front has evolved from the original distance excitation point position to anisotropy with axial symmetry to a regular circle, while the transverse wave no longer appears, and the slow longitudinal wave characteristics remain basically unchanged. This is consistent with the wave field characteristics of the two-phase isotropic medium. Combined with the characteristics of the energy changes before and after, it can be judged that on this basis, the anisotropy parameter is taken to zero, and the medium degenerates into a two-phase viscoelastic isotropic medium. Again, let the longitudinal and transverse wave quality factors Q P and Q Sis infinite, the wavefront shape in the wave field is consistent, but the energy is enhanced, the longitudinal and transverse wave quality factors are infinite, and the seismic wave energy attenuation phenomenon caused by the viscosity of the medium no longer occurs, which conforms to the wave field characteristics of a two-phase isotropic perfect elastic medium, and the medium degenerates into a two-phase isotropic perfect elastic medium. Finally, let the relevant parameters describing the two-phase medium be = = = = =0, the wave field energy is further enhanced, and the slow longitudinal wave no longer exists. This is because the parameter describing the interaction between the solid and liquid phases of the two-phase medium is zero and no longer plays a role in the wave equation. The forward modeling result of the medium no longer has the characteristics of a two-phase medium, but evolves into a uniform single-phase isotropic perfectly elastic medium.

[0208] Figure 3 is relative Figure 2 Numerical simulations of sequence-modified degradation were performed to validate the forward modeling results, and Figure 2 The forward modeling series uses the same set of parameter tables. The idea of ​​this forward modeling is to first degenerate the two-phase viscoelastic TTI medium into a single-phase viscoelastic TTI medium, and then further eliminate the anisotropy and viscoelasticity of the medium, and finally degenerate into a single-phase homogeneous isotropic perfect elastic medium. The specific steps are: from left to right, first assign the parameters describing the two-phase properties of the medium = = = = =0, the wave field snapshot has no slow longitudinal wave component, while other components are retained. The wave field morphology of fast longitudinal wave and transverse wave has not changed, and their symmetry is consistent, but the energy has been significantly enhanced. This is because there is no effect of the two-phase parameter, the medium does not show biphasicity, and there is no energy attenuation caused by the interaction and coupling effect between solid and fluid. Therefore, after this treatment, the medium is degenerated from a two-phase viscoelastic TTI medium to a single-phase viscoelastic TTI medium; secondly, the polarization angle θ=0, both fast longitudinal wave and transverse wave splitting exist, and the overall energy change of the two components is not obvious, but the wave front morphology rotates, and the symmetry axis changes from the original tilt of 45° to coincide with the central vertical line. Such a wave field conforms to the wave field characteristics of a single-phase VTI medium. Combined with the invariance of the medium energy, it can be judged that the medium has evolved from a single-phase viscoelastic TTI medium to a single-phase viscoelastic VTI medium; thirdly, based on the previous forward modeling, the anisotropy parameter = = = 0, the seismic wave field has only one type of wave, namely, longitudinal wave (since the forward simulation of this application example uses a longitudinal wave source, no transverse wave will be generated in a homogeneous isotropic medium without separating the longitudinal and transverse wave fields), and the wave field energy difference is not large, which shows that the medium degenerates into a single-phase viscoelastic medium; finally, the longitudinal and transverse wave quality factors Q P and Q S When the value is infinite, the medium degenerates into a uniform isotropic perfectly elastic medium, the wave field morphology does not change, but the wave field energy is significantly enhanced, which conforms to the wave field characteristics of an isotropic perfectly elastic medium. Figure 2 and Figure 3 From the results, it can be seen that as long as the values ​​of relevant parameters are adjusted, the forward modeling results of various corresponding types of media can be obtained based on the finite difference forward simulation of the present invention, and finally the wave field results of isotropic perfectly elastic media can be degenerated, thereby verifying the accuracy and comprehensiveness of the wave equation.

[0209] (b) Optimization of two-way wave illumination and observation system for the depression work area model

[0210] A two-phase viscoelastic TTI medium model is established based on a depression model of a certain work area. The model grid size is 288 151, the grid spacing is = =10 m, multi-shot excitation was carried out on the surface, the number of shots was 20, the shot spacing was 140 m, the source was uniformly Ricker wavelet longitudinal wave source, the main frequency was 25 Hz, and the time sampling interval of the forward simulation was 0.5 ms. The distribution of the parameter values ​​of the medium model is shown in Figure 4 , the first row from left to right is , , , the second row from left to right is , , , the third row from left to right is , , .like Figure 5A and Figure 5B As shown in the figure, the illumination results of the solid and liquid phases show abnormal values, among which the illumination energy of the liquid phase shows a negative number, which is obviously inconsistent with reality. Figure 5C The source illumination result calculated by the two-phase illumination operator is normal in shape. This result shows the correctness of the two-phase illumination operator, and the subsequent reverse illumination and optimized encryption results can further prove this point. Fig. 6A is the source lighting result, Figure 6B is the reverse illumination result of the model. Fig. 6AThe illumination result of the source direction is intercepted, and the illumination result of the area below 400 m is cut. The normalized isoline of illumination energy is drawn. The shadow area is located according to the isoline, and the excitation area and excitation position of the reverse illumination are determined. The isoline is shown in Figure 7 . Figure 8 The normalized curve of the surface energy of the reverse illumination is shown. From the figure, it can be determined that the energy maximum range is roughly around several wave peaks between 200 and 700 m. The energy in the horizontal distance between 1500 and 2800 m is very low. For energy optimization of the target area, it is possible to consider reducing the density of shot points in this range. According to the reverse illumination positioning encryption interval, the shot points are encrypted in the encrypted interval, and the number of shots in the low-energy interval indicated by the reverse illumination surface energy curve is correspondingly reduced. The encrypted source illumination results are shown in Fig. 9 . Fig.10 The comparison results of the normalized illumination energy of the shadow area of ​​the medium model before and after encryption are shown. From the reflection of this result, it can be seen that the illumination energy of the shadow area of ​​the observation system after the densification of the shot points has been significantly enhanced. This result conforms to the reciprocal principle. Therefore, the observation system is optimized, and the correctness of the two-phase illumination operator and the two-phase viscoelastic TTI medium wave equation are verified.

[0211] (c) Optimization of the two-way wave illumination and observation system of the Marmours model

[0212] The Marmours model of two-phase viscoelastic TTI properties is established, and the model grid size is 298 170, the grid spacing is = =10 m, multi-shot excitation was carried out on the surface, the number of shots was 20, the shot spacing was 150 m, the seismic source was uniformly Ricker wavelet longitudinal wave source, the main frequency was 25 Hz, and the time sampling interval of the forward simulation was 0.5 ms. The distribution of the parameter values ​​of the medium model is shown in Fig.11 , the first row from left to right is , , , the second row from left to right is , , , the third row from left to right is , , . Fig. 12A It is the result of source illumination, intercepting the area below 1200 m in depth, and locating the shadow area according to the distribution of regional energy. The energy value distribution is shown in Fig. 12B . Fig.13A is the reverse illumination result of the Marmours model, Fig. 13BThe normalized curve of reverse illumination surface energy is shown in the figure. The high energy interval is roughly around several wave peaks between 1250 and 2100 m, and the low energy interval is within the horizontal distance of 0 to 1100 m. The number of shots in the high energy interval of reverse illumination surface energy is increased, and the number of shots in the low energy interval is reduced. The source illumination result after the shot point encryption is shown in Fig.14A and 14B . Fig.15 The comparison results of the normalized illumination energy of the shadow area of ​​the Marmours model before and after encryption are shown. The illumination energy of the shadow area has been significantly enhanced after encryption. This result also conforms to the reciprocity principle, which shows that the observation system has been optimized. At the same time, the correctness and applicability of the illumination operator and wave equation have been verified.

[0213] The present invention combines the anisotropy theory of TI media with the Carcione viscoelastic anisotropy theory to describe the constitutive relationship of solid media, and applies its mechanism to the Biot two-phase medium system, and successfully derives the wave equation of the two-phase viscoelastic TTI medium. The wave equation is numerically simulated using the high-order staggered grid finite difference method, and the results clearly reflect the various wave field information of the two-phase viscoelastic TTI medium. The numerical simulation results include fast longitudinal waves, shear waves and slow longitudinal waves, and its wave field characteristics have the characteristics of TTI media, viscoelastic media and two-phase media at the same time, which is consistent with the relevant theoretical predictions. By comparing the forward modeling of seismic waves of various types of media, the wave field characteristics of the two-phase viscoelastic TTI medium can be clearly demonstrated, verifying the validity of the derived wave equation. The seismic record curve generated by the forward modeling is used to qualitatively study the influence of relevant parameters on the wave field of the two-phase viscoelastic TTI medium, thereby verifying the relevant theory to a certain extent. In addition, the two-way wave illumination and observation system optimization of the two-phase viscoelastic TTI medium are performed using the newly established two-phase illumination operator. The results show that the calculation method of the two-phase illumination operator can well verify the principle of seismic wave reciprocity and the observation system optimization principle of the statistical illumination method. The illumination analysis process demonstrates the application of the two-phase viscoelastic TTI medium wave equation in complex models, and further verifies the correctness and applicability of the wave equation.

[0214] Of course, the above description is not a limitation of the present invention, and the present invention is not limited to the above examples. Changes, modifications, additions or substitutions made by technicians in this technical field within the essential scope of the present invention should also fall within the protection scope of the present invention.

Claims

1. A seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium, characterized in that: The steps include: Step S1, introducing Biot two-phase medium theory to describe the two-phase medium; Step S2, using the TTI medium as a solid skeleton, deriving the constitutive equation of the TTI medium; Step S3, introducing Carcione's viscoelastic anisotropic medium model and combining it with the Boltzmann superposition principle to obtain the stress-strain relationship of the solid skeleton; Step S4, embedding the viscoelastic constitutive relation of the solid phase into Biot's two-phase medium theory, and deriving a wave equation applicable to the two-phase viscoelastic TTI medium; Step S5, using a high-order staggered grid finite difference method to perform seismic wave forward simulation on the derived first-order velocity stress equation of the two-phase viscoelastic TTI medium; Step S6, establishing an illumination operator, defining source direction and reverse illumination operators, calculating source direction illumination through forward simulation, and analyzing the illumination effects on the target area under different configurations; Step S7: Based on the reciprocity principle, analyze the reverse illumination results, identify the surface locations that contribute more to the target area, and optimize the excitation positions and shot point density to achieve optimal coverage.

2. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 1, characterized in that: In step S2, the constitutive equation of the TTI medium is: ； In the formula, is the TTI medium stress vector, is the TTI medium strain vector.

3. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 2, characterized in that: In step S3, according to Carcione viscoelastic anisotropy theory and Boltzmann superposition principle, the stress-strain relationship of the solid skeleton is: ； In the formula, is the stress vector, is the first-order time derivative of the strain vector, represents the relaxation matrix.

4. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 3, characterized in that: In step S4, the wave equation applicable to the two-phase viscoelastic TTI medium is: ； In the formula, is the displacement vector; express The second derivative with respect to time; is the body force vector, where , , represents the solid phase fraction, , , represents the liquid phase fraction, Compared to A 7×7 matrix with a row vector and a column vector with 0 elements is added respectively; is the density matrix; is a partial differential operator; is the dissipation coefficient matrix; is the sum and difference operator.

5. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 4, characterized in that: In step S5, the first-order time derivatives of the stress and strain of the wave equation of the two-phase viscoelastic TTI medium are obtained to obtain the first-order velocity stress equation of the two-phase viscoelastic TTI medium; The first-order time derivative equation of the memory variable is: ； In the formula, relaxation time Calculated by quality factor, the calculation formula is: ； In the formula, , and are the longitudinal wave quality factor and the shear wave quality factor, It is the reciprocal of the main frequency of the source wavelet.

6. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 5, characterized in that: Taking into account the viscoelasticity and TI anisotropy of the medium, the necessary stability condition for the forward simulation of the first-order velocity stress equation of the two-phase viscoelastic TTI medium using the high-order staggered grid finite difference method is: ； In the formula, .

7. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 6, characterized in that: In step S6, the illumination operator expression is: Source lighting operator: ； In the formula, is the excitation position of the observation system, In the two-phase viscoelastic medium model region Coordinates of any point; The reverse illumination operator is: ； In the formula, are the coordinates of each position in the observation range, The source illuminates the target point in the shadow area. The shaded area.

8. The seismic wave forward modeling and two-way wave illumination method for a two-phase viscoelastic TTI medium according to claim 7, characterized in that: In step S7, the illumination function of the optimized and encrypted observation system for the two-phase medium is: ； In the formula, It is an encrypted surface observation system.

Images