A method for simulating seismic wave propagation characteristics of porous rock based on multi-mechanism coupling

Abstract

The application belongs to the technical field of seismic exploration, and is a pore rock seismic wave propagation characteristic simulation method based on multi-mechanism coupling, comprising: simulating the propagation of seismic waves in viscoelastic porous rocks by using a seismic wave propagation model; wherein the modeling process of the seismic wave propagation model comprises: coupling the squirt flow mechanism and the Biot theory framework; using the Zener model to represent the inherent viscoelasticity of the solid skeleton, and coupling the squirt flow mechanism and the Biot theory framework to obtain the total bulk modulus and the total shear modulus; calculating the memory variable according to the total bulk modulus and the total shear modulus; obtaining the time-domain relaxation function containing the memory variable; substituting the time-domain relaxation function containing the memory variable into the frequency-domain stress-strain relationship to obtain the time-domain velocity-stress equation set containing all memory variables; and combining the first-order differential equation of the memory variable to obtain the seismic wave propagation model. The application can accurately and comprehensively simulate the wave field characteristics such as reflection and diffraction in complex heterogeneous media.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of seismic exploration, and in particular to a simulation method for seismic wave propagation characteristics of porous rock based on multi-mechanism coupling. BACKGROUND

[0002] In the field of geophysical exploration, the propagation characteristics of seismic waves in porous and dissipative media are the key basis for understanding the underground structure and resource distribution. When seismic waves propagate in fluid-saturated rocks and other porous media, they exhibit significant viscoelastic characteristics, specifically embodied in frequency-dependent energy loss (attenuation) and velocity dispersion, which are of great significance to environmental monitoring, geothermal exploration, oil and gas exploitation, geotechnical engineering, and seismic research, and many other applications of earth science.

[0003] Fluid-saturated rock is usually modeled as a "solid skeleton-pore fluid" two-phase system, where the fluid fills the pore space of the solid skeleton. The classical theory proposed by Biot first systematically describes the propagation of seismic waves in such a two-phase system, and clearly identifies the flow of the fluid relative to the solid skeleton as one of the key mechanisms leading to seismic wave attenuation. This theory has since become the basis for understanding the wave propagation characteristics in sedimentary rocks.

[0004] However, existing research shows that seismic wave attenuation and velocity dispersion are significantly present in a wide frequency range, and the assumption of uniform porosity in Biot's theory is difficult to capture the mechanical response of the complex pore structure of real rocks. Wave-induced fluid flow, as a key mechanism of viscoelastic attenuation in fluid-saturated rock, plays a role on multiple scales, with microscale jet flow being particularly important. This mechanism is due to fluid exchange between pores of different compressibility or permeability, and although compliant pores only account for a small fraction of the total pore space, they have a significant impact on the overall mechanical properties of the rock. For the jet flow mechanism, the academic community has developed a simplified and efficient analytical model by introducing a jet flow aspect ratio parameter, which is used to evaluate the seismic dispersion and attenuation characteristics of isotropic porous rocks. However, the fluid-related dissipation mechanism cannot fully explain the energy loss in real rocks, and the intrinsic viscoelastic properties of the rock solid skeleton itself are not ideal elastic bodies, which also have a significant impact on seismic wave attenuation and dispersion. However, existing models have not fully integrated this key mechanism. SUMMARY

[0005] The embodiments of the present application provide a simulation method for seismic wave propagation characteristics of porous rock based on multi-mechanism coupling, which solves the problem that the coupling effect of three key attenuation mechanisms in a wide frequency range cannot be captured in a unified framework, resulting in incomplete simulation of attenuation peak and dispersion characteristics, and difficulty in meeting the needs of complex geological scenarios.

[0006] According to the simulation method for seismic wave propagation characteristics of porous rock based on multi-mechanism coupling provided by the embodiments of the present application, the method comprises:

[0007] constructing a seismic wave propagation model;

[0008] simulate propagation of seismic waves in viscoelastic porous rock using the seismic wave propagation model;

[0009] The process of constructing the seismic wave propagation model comprises:

[0010] obtain a stress-strain relationship in frequency domain by coupling the jet-flow mechanism with the Biot theory framework;

[0011] convert the frequency-dependent complex modulus in the stress-strain relationship in frequency domain into a complex modulus GZM model representation form; the frequency-dependent complex modulus comprises a bulk modulus corresponding to a mid-frequency attenuation characteristic of fluid flow between cracks and a shear modulus corresponding to a mid-frequency attenuation characteristic of fluid flow between pores; the complex modulus GZM model representation form comprises a GZM model representation form of the bulk modulus and a GZM model representation form of the shear modulus;

[0012] obtain a total bulk modulus and a total shear modulus, and couple the total bulk modulus and the total shear modulus with the GZM model representation form of the bulk modulus and the GZM model representation form of the shear modulus respectively to represent inherent viscoelasticity of the solid skeleton by using a Zener model, so as to obtain the total bulk modulus and the total shear modulus;

[0013] calculate a memory variable according to the total bulk modulus and the total shear modulus;

[0014] multiply a high-frequency unrelaxed modulus by an instantaneous strain rate, and subtract a sum of all memory variables to obtain a time-domain relaxation function containing the memory variable;

[0015] substitute the time-domain relaxation function containing the memory variable into the stress-strain relationship in frequency domain, replace a product of the complex modulus and the strain in the frequency domain with a product of the high-frequency unrelaxed modulus and the instantaneous strain rate in the time domain, and subtract the sum of all memory variables to obtain a time-domain stress-strain relationship containing the memory variable, and then derive the time-domain stress-strain relationship with respect to time to obtain a time-domain velocity-stress equation group containing all memory variables; combine a first-order differential equation of the memory variable to obtain the seismic wave propagation model.

[0016] Further, the stress-strain relationship in frequency domain obtained by coupling the jet-flow mechanism with the Biot theory framework comprises:

[0017] ;

[0018] ;

[0019] wherein, denotes a total stress, denotes a total fluid pressure; denotes the solid strain tensor, is the Kronecker symbol, when , = 1, otherwise 0; denotes the total fluid pressure, denotes the solid phase volumetric strain, denotes the fluid phase relative displacement divergence, is the representation of the shear modulus as a function of the angular frequency, the frequency-dependent behavior of G represents the strength of the correlation between the fluid pressure and the deformation at different frequencies, is the coupling modulus, is the effective stress coefficient.

[0020] Further, the bulk modulus is calculated from the effective fluid modulus, and the shear modulus is calculated from the bulk modulus, and the formula is represented as:

[0021] ,

[0022] ,

[0023] wherein, is the bulk modulus, is the bulk modulus of the solid, represents the bulk modulus of the fluid, is the shear modulus of the solid, represents the porosity; is the bulk modulus of the solid particles, is the shear modulus, is the effective fluid modulus.

[0024] Further, the formula for representing the effective fluid modulus is: wherein, represents the shear viscosity of the fluid, is the jet flow aspect ratio is the parameter related to the pressure diffusion of the jet flow, , is the effective fluid modulus, is the imaginary unit, is the angular frequency, is the shear viscosity of the fluid, denotes the bulk modulus of the fluid, is the complex modulus of the jet flow mechanism.

[0025] Further, the frequency-dependent complex modulus in the frequency domain stress-strain relationship is converted into the complex modulus GZM model representation, including:

[0026] The complex modulus of the jet flow mechanism is approximated by a rational function, and the optimal coefficients of the rational function are obtained by genetic algorithm, so that the error between the optimal coefficients and the theoretical values in the dominant frequency band of the jet flow is small.

[0027] The effective fluid modulus is calculated by the rational function with the optimal coefficients, and the bulk modulus and shear modulus are calculated by the effective fluid modulus.

[0028] The parameters of the complex modulus of the GZM model are optimized by genetic algorithm, so that the error between the complex modulus of the GZM model and the bulk modulus or shear modulus is minimized, and the parameters of the complex modulus of the GZM model corresponding to the bulk modulus or shear modulus are obtained.

[0029] Based on the optimized parameters of the complex modulus of the GZM model, the GZM model expression of the bulk modulus and the GZM model expression of the shear modulus are obtained.

[0030] Further, the Zener model is used to represent the inherent viscoelasticity of the solid skeleton, which is represented as:

[0031] = , wherein represents the skeleton strain relaxation time, is the skeleton stress relaxation time, is the angular frequency, is the imaginary unit, is the Zener model factor of the inherent viscoelasticity of the solid skeleton;

[0032] ;

[0033] ;

[0034] is the skeleton mass factor, is the center frequency.

[0035] Further, the total bulk modulus is the product of the Zener model factor of the inherent viscoelasticity of the solid skeleton and the GZM model expression of the bulk modulus;

[0036] The total shear modulus is the product of the Zener model factor of the inherent viscoelasticity of the solid skeleton and the GZM model expression of the shear modulus.

[0037] Further, the memory variable is represented as: , is the time variable, is the deformation rate function varying with time, is the relaxation strength of the first variable, is the time derivative of the first variable, =3, For the first The stress relaxation time of each variable.

[0038] Furthermore, the first-order differential equation of the memory variable is: The current memory variable is updated by comparing the current strain with the memory variable value from the previous time step. Based on the first-order differential equation, the values ​​of all memory variables are updated at each time step. The updated memory variables are then substituted into the time-domain velocity-stress equation set containing all memory variables to calculate the stress and fluid pressure at the current time step. The calculated stress and fluid pressure are then used to update the velocity and displacement fields of the solid and fluid to obtain the strain field at the next time step. This process is repeated until the simulation is complete.

[0039] Compared with the prior art, the advantages of this application are as follows: the method of this application can accurately and comprehensively simulate the wave field characteristics such as reflection and diffraction in complex heterogeneous media, and the multi-mechanism attenuation effect is basically consistent with the actual strata, meeting the needs of complex geological scenarios. Attached Figure Description

[0040] Figure 1 A flowchart illustrating the method provided in the embodiments of this application;

[0041] Figure 2 Provided for the embodiments of this application Wavefield snapshots and seismic records, where (a) represents the solid particle velocity components in the Marmousi model. (a) Wavefield snapshots and seismic records, (b) showing the solid particle velocity components in the Marmousi model. Wavefield snapshots and seismic records;

[0042] Figure 3 Provided for the embodiments of this application Wavefield snapshots and seismic records, where (a) represents the solid particle velocity components in the Marmousi model. (a) Wavefield snapshots and seismic records, (b) showing the solid particle velocity components in the Marmousi model. Wavefield snapshots and seismic records. Detailed Implementation

[0043] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the scope of this application.

[0044] See Figure 1 As shown in the embodiment of this application, a method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling is provided, including: constructing a seismic wave propagation model;

[0045] A seismic wave propagation model was used to simulate the propagation of seismic waves in viscoelastic porous rocks;

[0046] The modeling process of the seismic wave propagation model includes:

[0047] By coupling the jet flow mechanism with the Biot theoretical framework, the frequency domain stress-strain relationship was obtained;

[0048] The frequency-dependent complex modulus in the frequency domain stress-strain relationship is transformed into a GZM model representation of the complex modulus; the frequency-dependent complex modulus includes the bulk modulus and the shear modulus, the bulk modulus corresponds to the mid-frequency attenuation characteristics of fluid flow in the fracture, and the shear modulus corresponds to the mid-frequency attenuation characteristics of fluid flow in the pore; the GZM model representation of the complex modulus includes the GZM model representation of the bulk modulus and the GZM model representation of the shear modulus.

[0049] The Zener model is used to characterize the intrinsic viscoelasticity of the solid skeleton. It is coupled with the GZM model representations of the bulk modulus and the shear modulus to obtain the total bulk modulus and the total shear modulus.

[0050] Calculate the memory variables based on the total bulk modulus and total shear modulus;

[0051] The product of the high-frequency unrelaxed modulus and the instantaneous strain rate is then subtracted from the sum of all memory variables to obtain the time-domain relaxation function containing memory variables.

[0052] Substituting the time-domain relaxation function containing memory variables into the frequency-domain stress-strain relationship, the product of the complex modulus and strain in the frequency domain is replaced in the time domain with the product of the high-frequency unrelaxed modulus and the instantaneous strain rate. Subtracting the sum of all memory variables yields the time-domain stress-strain relationship containing memory variables. Differentiating the time-domain stress-strain relationship with respect to time yields the time-domain velocity-stress equations containing memory variables. Combining these with the first-order differential equations of the memory variables, the seismic wave propagation model is obtained.

[0053] The core of the jet flow mechanism is the wave-induced fluid flow between microscopic pores. Its energy dissipation and dispersion characteristics are closely related to pore geometry and frequency. The influence of pore geometry is parameterized by introducing the jet flow aspect ratio, where the jet flow aspect ratio is: , The pore size is the fracture diameter. For the fluid pressure diffusion path length associated with the jet flow, when When the pore size is smaller, the pore geometry is closer to a narrow, elongated fissure, resulting in a slower fluid pressure diffusion rate and thus enhancing the dispersion effect. The aspect ratio of the jet flow is used to describe a one-dimensional linear diffusion path, which better matches the topology of real rock fractures.

[0054] The jet flow causes fluid to flow between cracks and pores. Its impact on rock stiffness needs to be quantified using the effective fluid modulus, combined with the fluid's shear viscosity and angular frequency. The effective fluid modulus is derived from the aspect ratio of the jet flow: ,in, Represents the shear viscosity of a fluid. The aspect ratio of the jet flow Parameters related to jet pressure diffusion, , For effective fluid modulus, The imaginary unit, Angular frequency, For fluid shear viscosity, It represents the bulk modulus of a fluid.

[0055] Substituting the obtained jet flow aspect ratio into the above formula yields the effective fluid modulus, which includes terms related to the jet flow aspect ratio. This directly affects the frequency dependence. As the frequency increases, the imaginary part of the effective fluid modulus increases, and the energy dissipation is enhanced, which is consistent with the attenuation characteristics of the jet flow in the mid-frequency band. This allows for the simulation of the propagation speed and attenuation characteristics of seismic waves.

[0056] Based on the effective fluid modulus obtained from the above steps, the frequency-dependent bulk modulus and shear modulus in fluid-saturated rock are calculated using the following formulas:

[0057] ,

[0058] ,

[0059] Bulk modulus The bulk modulus of a solid. The bulk modulus of a fluid. For the shear modulus of a solid, Represents porosity; Bulk modulus of solid particles It is the shear modulus. This is the effective fluid modulus.

[0060] Based on the Biot theory framework and combining the aforementioned frequency-varying modulus (referring to the bulk modulus and shear modulus related to frequency variation), a frequency-domain stress-strain relationship integrating the jet flow mechanism is constructed. This relationship couples the jet flow mechanism with the Biot theory framework in the frequency domain. The total stress is simultaneously affected by solid elastic deformation, fluid pressure, and jet dissipation. The frequency-domain stress-strain relationship is as follows:

[0061] ;

[0062] ;

[0063] in, Indicates the total stress. Indicates the total fluid pressure; Represents the strain tensor of a solid. For Kroneck symbol, when hour, =1, otherwise 0; Indicates the total fluid pressure. Represents the volumetric strain of the solid phase. Represents the divergence of relative fluid displacement. It is a representation of the shear modulus as a function of angular frequency. The frequency-varying characteristics represent the strength of the correlation between fluid pressure and deformation at different frequencies. For coupling modulus, The effective stress coefficient;

[0064] In the frequency domain stress-strain relationship, there exists a variation with angular frequency. The changing physical parameters, their values ​​and characteristics, change with the frequency of seismic waves. These parameters are directly related to the frequency dependence of the frequency-domain stress-strain relationship, thus determining the attenuation and velocity dispersion characteristics of seismic waves, specifically including:

[0065] ,

[0066] ,

[0067] ,

[0068] ,

[0069] in, It is one of the basic parameters describing the elastic deformation of a solid skeleton, reflecting the relationship between stress and strain under the jet flow mechanism; Porosity;

[0070] The coupling modulus represents the strength of the coupling between the fluid and the solid skeleton, quantifying the interaction between fluid pressure and solid deformation. As frequency changes, It also exhibits frequency-varying characteristics, reflecting the differences in the efficiency of fluid-solid interaction at different frequencies;

[0071] The effective stress coefficient represents the proportion of pore pressure to the total stress. It is a key parameter connecting the total stress and pore pressure. At low frequencies, the fluid can flow freely, and the pore pressure is transmitted more fully. When the value is large and the frequency is high, fluid flow is restricted, and the effect of pore pressure is weakened. Follow Increase and decrease;

[0072] The frequency-dependent characteristic represents the modulus of the correlation strength between fluid pressure and deformation at different frequencies;

[0073] It is a representation of the change of bulk modulus with angular frequency. It is a representation of the shear modulus as a function of angular frequency.

[0074] By integrating the Biot global flow and jet flow mechanisms in the above manner, two significant attenuation peaks are observed in the wide frequency band, corresponding to different physical mechanisms: the jet flow peak in the mid-frequency range and the Biot peak in the high-frequency range. Compared with the traditional model, the integrated jet flow peak is located lower, has a higher degree of agreement with the numerical simulation results, and is closer to the real fracture-pore geometry.

[0075] The above process is in the frequency domain, but actual seismic wave propagation simulation needs to be carried out in the time domain and needs to be applied to time domain simulation.

[0076] However, traditional time-domain implementation methods have the problem that the frequency-varying modulus of the jet flow mechanism contains a complex hyperbolic function, which cannot be directly converted into a time-domain expression. Furthermore, time-domain convolution operations have high computational cost and memory requirements, making it difficult to efficiently achieve large-scale simulations.

[0077] To address this issue, a generalized Zener model (GZM model) combined with a genetic algorithm is used to approximate the frequency-dependent complex modulus corresponding to the jet flow mechanism. This approximates the frequency-dependent complex modulus corresponding to the jet flow mechanism (i.e.,...) , and Transforming the data into a memory variable form suitable for time-domain simulation, specifically including:

[0078] The GZM model describes the viscoelastic properties of materials by superimposing multiple Zener elements (relaxation mechanisms). The complex modulus of the GZM model... for: ,in, For low-frequency relaxation modulus, For the number of Zener components, and The first The strain relaxation time and stress relaxation time of a Zener element; the Zener element is the basic unit for describing the viscoelastic relaxation mechanism of a single frequency band.

[0079] To address the difficulty of time-domain transformation of the frequency-varying modulus of the jet flow mechanism, this paper transforms the frequency-dependent complex modulus in the frequency-domain stress-strain relationship into a complex modulus representation using the GZM model. This includes: approximating the complex modulus of the jet flow mechanism using rational functions through a genetic algorithm. (terms), in rational function form, are defined as: ,use Approximately obtained , , , , The coefficients obtained by the genetic algorithm optimization are the optimal coefficients with the smallest error compared to the theoretical values ​​in the dominant jet frequency band after executing the genetic algorithm. , , , Substitute into rational functions ;

[0080] The effective fluid modulus is calculated using the rational function corresponding to the optimal coefficients. The bulk modulus and shear modulus are then calculated using the effective fluid modulus, including: substituting the optimized rational function into the expression for the effective fluid modulus. The bulk modulus and shear modulus of the jet were obtained.

[0081] The parameters of the complex modulus of the GZM model are optimized using a genetic algorithm to minimize the error between the complex modulus and the bulk modulus or shear modulus of the GZM model. Specifically, the parameters of the complex modulus of the GZM model are obtained by optimizing the complex modulus of the GZM model using a genetic algorithm. middle , , Parameters such as these make the complex modulus of the GZM model... The error with the frequency-varying modulus of the jet stream is minimal. The frequency-varying modulus includes bulk modulus and shear modulus.

[0082] Parameters of the complex modulus based on the optimized GZM model ( , , By utilizing its superposition property of multiple relaxation mechanisms, the stress-strain relationship containing convolution operations in the time domain is transformed into a set of differential equations about memory variables, thereby achieving efficient simulation of jet flow effects in the time domain;

[0083] The time-domain relaxation function of the GZM model is : ,in, For unit transition function, For high-frequency unrelaxed modulus, For the first The relaxation strength of each Zener element. To avoid temporal convolution operations, the temporal relaxation function of the GZM model is transformed into an equation with memory variables to represent the stress relaxation time. In this case, the number of Zener elements is two. Two Zener elements can accurately match the attenuation and dispersion characteristics of the frequency-varying modulus in the dominant jet frequency band. Too many Zener elements would increase the computational load and memory requirements of the time-domain simulation. Two Zener elements balance accuracy and efficiency, transforming the temporal relaxation function of the GZM model... The equation can be transformed into an equation containing memory variables as follows:

[0084] = = ,

[0085] It is a deformation rate function that varies with time. It is the time derivative of the memory variable.

[0086] Substituting the time-domain relaxation functions, converted to memory variable form, into the frequency-domain stress-strain relationship yields a time-domain stress-strain relationship containing memory variables that can be directly used for time-domain simulation. However, this is the case with two Zener elements.

[0087] For each frequency-varying modulus characterized by the Zener model, based on the corresponding stress relaxation time... and relaxation strength Establish a memory variable for each Zener element. .

[0088] The core of transforming the frequency-domain stress-strain relationship into a time-domain form lies in replacing the product of the complex modulus and complex strain in the frequency domain with the product of the high-frequency unrelaxed modulus and the instantaneous strain rate in the time domain. Then, subtracting the sum of all relevant memory variables yields a time-domain stress-strain relationship including memory variables. Differentiating this time-domain stress-strain relationship with time variables further yields a time-domain velocity-stress equation system including memory variables. The maximum value is 2:

[0089] ,

[0090] ,

[0091] ,

[0092] The formulas in the first and second rows are the time derivatives of the total stress tensor. , , , , , The first line represents the time derivatives of the different planar components of the stress tensor, describing the rate of change of stress over time. It is a core output of the finite-difference time-domain simulation, directly reflecting the energy dissipation during wave propagation; the second line represents the total fluid pressure. time derivative ; These are the velocity vector components of the solid particles; It is the solid velocity relative to spatial coordinates The partial derivatives; It is the volumetric strain rate of the solid; It is the relative volumetric strain rate of the fluid; It is the plane shear strain rate of ij. for , as well as , for , as well as ,

[0093] It is the high-frequency unrelaxed modulus of the bulk modulus of the jet flow GZM model; It is the high-frequency unrelaxed modulus of the shear modulus in the jet flow GZM model; It is the effective stress coefficient The high-frequency unrelaxed modulus; yes The high-frequency unrelaxed modulus; It is the coupling modulus The high-frequency unrelaxed modulus;

[0094] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the planar direction;

[0095] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the planar direction;

[0096] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the planar direction;

[0097] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the direction;

[0098] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the direction;

[0099] It is the first to integrate Biot's global flow and jet flow mechanisms. Each Zener element corresponds to The memory variable for the shear modulus in the direction;

[0100] It is the first Each Zener element corresponds to a coupling modulus Memory variables;

[0101] It is the first to integrate Biot's global flow and jet flow mechanisms. The memory variable for the bulk modulus of a Zener element;

[0102] It is the first Effective stress coefficient of each Zener element Memory variables;

[0103] It is the first Zener element The memory variable of the modulus.

[0104] The Zener model is used to characterize the intrinsic viscoelasticity of the solid skeleton. It is coupled with the GZM model representations of the bulk modulus and shear modulus to obtain the total bulk modulus and total shear modulus, including:

[0105] The Zener model is introduced to characterize the intrinsic viscoelasticity of the solid framework. A new memory variable is added to achieve coupling with the Biot framework and the jet flow mechanism through three mechanisms, as detailed below:

[0106] The Zener model is used to characterize the intrinsic viscoelasticity of the solid skeleton, incorporating the skeleton mass factor and the skeleton center frequency. ,in, Indicates the strain relaxation time of the skeleton. The skeleton stress relaxation time; the skeleton quality factor is The center frequency is The core parameters for calculating the Zener model, which is used to characterize the intrinsic viscoelasticity of a solid skeleton, are the skeleton strain relaxation time and the skeleton stress relaxation time.

[0107] Based on the jet flow mechanism, the bulk modulus and shear modulus are used to obtain the GZM model expression after optimization by a genetic algorithm. and The fracture jetting effect is integrated into the skeleton model, retaining two relaxation elements of the jetting mechanism, corresponding to the mid-frequency attenuation characteristics of fluid flow between fractures and pores. The total bulk modulus and total shear modulus, including multiple mechanisms, are calculated by the following formula, which is a direct result of this coupling action. That is, the coupling of the jetting flow and the skeleton model is completed first, and then the coupled total bulk modulus and total shear modulus are output: The bulk modulus refers to the overall viscoelasticity of the integrated jet stream and the inherent viscoelasticity of the solid skeleton. The total shear modulus refers to the sum of the intrinsic viscoelasticity of the jet stream and the solid skeleton. These are the bulk modulus and shear modulus that vary frequently due to the jet flow mechanism. The Zener model factor refers to the intrinsic viscoelasticity of a solid skeleton. , These refer to the high-frequency unrelaxed modulus of the bulk modulus and the high-frequency unrelaxed modulus of the shear modulus of the jet flow GZM model, respectively. , These refer to the relaxation strengths of the bulk modulus and shear modulus of the jet flow GZM model, respectively. The contribution weight of each Zener element to volume and shear attenuation is as follows: the larger the value, the stronger the attenuation in that frequency band. Strain relaxation time of bulk modulus in the jet flow GZM model It is the strain relaxation time of the shear modulus in the jet flow GZM model; Indicates the strain relaxation time of the skeleton. Given the skeleton stress relaxation time, the total bulk modulus and total shear modulus are expressed as follows:

[0108] ,

[0109] .

[0110] To adapt to time-domain simulations, the coupled total modulus (total bulk modulus and total shear modulus) is transformed from a two-element Zener element to a GZM model expression containing three Zener elements (the first two elements come from the jet GZM model, and the third comes from the skeleton Zener model). The parameters of each Zener element are determined, including the high-frequency non-relaxation modulus (high-frequency limit value), stress and shear modulus relaxation time (characteristic time of each mechanism), and stress and shear modulus relaxation strength (contribution weight of each mechanism to the total attenuation).

[0111] The time-domain relaxation function of the GZM model after integrating the jet flow and the intrinsic viscoelasticity of the solid skeleton. The time-domain relaxation function representation containing memory variables is as follows: = = ,

[0112] Substituting the time-domain relaxation function, which integrates the intrinsic viscoelasticity of the jet flow and the solid skeleton and is converted into a memory variable form, into the frequency-domain stress-strain relationship, we obtain a time-domain stress-strain relationship containing all memory variables that can be directly used for time-domain simulation. However, this is the case with three Zener elements.

[0113] The first-order differential equation for the memory variable is established as follows: By updating the current memory variable with the value of the memory variable from the previous moment, the time-domain discretization of the stress-strain relationship is achieved. For time variables The changing rate of deformation function It is a time-varying variable Changing memory variables.

[0114] The memory variable is represented as: , It is a time variable. It is a time-varying variable The changing rate of deformation function For the first The relaxation strength of each variable, It is the first Time derivatives of each variable =3, For the first The stress relaxation time of each variable.

[0115] At this point, in the time-domain velocity-stress equations of the memory variable... The maximum value is 3, and the time-domain velocity-stress equations, which include all memory variables, are expressed in the following form:

[0116] ,

[0117] ,

[0118] ,

[0119] in, It is the high-frequency unrelaxed modulus of the total bulk modulus; It is the high-frequency unrelaxed modulus of the total shear modulus;

[0120] It is the first step after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the planar direction;

[0121] It is the first step after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the planar direction;

[0122] It is the first after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the planar direction;

[0123] It is the first after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the direction;

[0124] It is the first after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the direction;

[0125] It is the first after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. Each Zener element corresponds to The memory variable for the total shear modulus in the direction;

[0126] It is the first after integrating the inherent viscoelasticity of the jet stream and the solid skeleton. The memory variable for the total bulk modulus of each Zener element;

[0127] The current memory variable is updated by comparing the current strain with the memory variable value from the previous time step. Based on the first-order differential equation of the memory variable, the values ​​of all memory variables are updated at each time step. The updated memory variable is then substituted into the time-domain velocity-stress equation system containing the memory variable to calculate the stress and fluid pressure at the current time step. The calculated stress and fluid pressure are then used to update the velocity and displacement fields of the solid and fluid to obtain the strain field at the next time step. This process is repeated until the simulation is complete.

[0128] By transforming the convolution operation in the stress-strain relationship, which depends on the entire history, into a recursive operation that depends only on the state of the previous time step, efficient time-domain discretization of the stress-strain relationship is achieved.

[0129] To verify the wavefield characteristics in complex inhomogeneous media, the Marmousi model was used for simulation, such as... Figure 2 yes Wavefield snapshots and seismic records, among whichFigure 2 (a) represents the velocity component of the solid particles in the Marmousi model. Wavefield snapshots and seismic records, Figure 2 (b) represents the velocity components of solid particles in the Marmousi model. Wavefield snapshots and seismic records, Figure 3 yes Wavefield snapshots and seismic records, among which Figure 3 (a) represents the velocity component of the solid particles in the Marmousi model. Wavefield snapshots and seismic records, Figure 3 (b) represents the velocity components of solid particles in the Marmousi model. Wavefield snapshots and seismic records. The method provided in this application can accurately capture the reflection, transmission, and diffraction phenomena of P-waves and S-waves, forming complex interferometric wavefields at fault and stratigraphic interfaces, verifying the model's ability to capture complex tectonic wavefields; P-waves, S-waves, and multiple sets of reflected waves can be identified in seismic records, with no obvious numerical oscillations in the waveforms, verifying the stability of the finite difference time-domain algorithm in heterogeneous media, compared with... Figure 2 and Figure 3 The records show that, with the increase of time, the amplitude of the reflected wave in the later stage conforms to the attenuation law of multiple mechanisms with the propagation distance; indicating that the method of this application can accurately simulate the wave field characteristics such as reflection and diffraction in complex heterogeneous media, and the multi-mechanism attenuation effect is basically consistent with the actual strata.

[0130] The above description is merely a preferred embodiment of this application and is not intended to limit this application. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling, characterized in that, The method includes: Construct a seismic wave propagation model; A seismic wave propagation model was used to simulate the propagation of seismic waves in viscoelastic porous rocks; The process of constructing a seismic wave propagation model includes: By coupling the jet flow mechanism with the Biot theoretical framework, the frequency domain stress-strain relationship is obtained, including: ； ； in, Indicates the total stress. Indicates the total fluid pressure; Represents the strain tensor of a solid. For Kroneck symbol, when hour, =1, otherwise 0; Represents the volumetric strain of the solid phase. Represents the divergence of relative fluid displacement. It is a representation of the shear modulus as a function of angular frequency. The frequency-varying characteristics represent the strength of the correlation between fluid pressure and deformation at different frequencies. For coupling modulus, The effective stress coefficient, It is one of the basic parameters describing the elastic deformation of a solid skeleton, reflecting the relationship between stress and strain under the jet flow mechanism; The frequency-dependent complex moduli in the frequency domain stress-strain relationship are transformed into a GZM model representation. The frequency-dependent complex moduli include bulk modulus and shear modulus. The bulk modulus corresponds to the mid-frequency attenuation characteristics of fluid flow between fractures, and the shear modulus corresponds to the mid-frequency attenuation characteristics of fluid flow between pores. The GZM model representation of the complex modulus includes GZM model representations of both bulk modulus and shear modulus. The GZM model describes the viscoelastic properties of the material by superimposing multiple Zener elements. The complex modulus of the GZM model... for: ,in, For low-frequency relaxation modulus, For the number of Zener components, and The first The strain relaxation time and stress relaxation time of a Zener element; the Zener element is the basic unit for describing the viscoelastic relaxation mechanism of a single frequency band. The total bulk modulus and total shear modulus are obtained by characterizing the intrinsic viscoelasticity of the solid skeleton using the Zener model, and coupling it with the GZM model representations of the bulk modulus and the shear modulus, respectively, to obtain the total bulk modulus and total shear modulus. Calculate the memory variables based on the total bulk modulus and total shear modulus; The product of the high-frequency unrelaxed modulus and the instantaneous strain rate is then subtracted from the sum of all memory variables, including: the memory variables of the shear modulus corresponding to the three planar directions and three axial directions of the Zener element that integrates the Biot global flow and jet flow mechanisms, and the coupling modulus corresponding to the Zener element. Memory variables, memory variables of the bulk modulus of Zener elements integrating Biot's global flow and jet flow mechanisms, and effective stress coefficients of Zener elements. Memory variables and Zener elements The memory variable of the modulus is used to obtain the time-domain relaxation function containing the memory variable after integrating the intrinsic viscoelasticity of the jet flow and the solid skeleton; Substituting the time-domain relaxation function containing memory variables into the frequency-domain stress-strain relationship, the product of the complex modulus and strain in the frequency domain is replaced in the time domain with the product of the high-frequency unrelaxed modulus and the instantaneous strain rate. Subtracting the sum of all memory variables yields the time-domain stress-strain relationship containing memory variables. Differentiating the time-domain stress-strain relationship with respect to time yields the time-domain velocity-stress equations containing all memory variables. These memory variables include: the memory variable representing the total shear modulus of the Zener element in the three planar directions and three axial directions after integrating the intrinsic viscoelasticity of the jet stream and the solid skeleton; the memory variable representing the total bulk modulus of the Zener element after integrating the intrinsic viscoelasticity of the jet stream and the solid skeleton; and the memory variable representing the coupling modulus of the Zener element. Memory variables, effective stress coefficient of Zener elements Memory variables and Zener elements The memory variable of the modulus; combined with the first-order differential equation of the memory variable, the seismic wave propagation model is obtained.

2. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 1, characterized in that, The bulk modulus is calculated based on the effective fluid modulus, and the shear modulus is calculated based on the bulk modulus, expressed by the following formula: ， ， in, Bulk modulus The bulk modulus of a solid. The bulk modulus of a fluid. For the shear modulus of a solid, Represents porosity; Bulk modulus of solid particles It is the shear modulus. This is the effective fluid modulus.

3. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 2, characterized in that, The formula for expressing the effective fluid modulus is: ,in, Represents the shear viscosity of a fluid. The aspect ratio of the jet flow Parameters related to jet pressure diffusion, , The imaginary unit, Angular frequency, For fluid shear viscosity, This is the complex modulus of the jet flow mechanism.

4. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 3, characterized in that, Transforming the frequency-dependent complex modulus in the frequency domain stress-strain relationship into a GZM model representation of the complex modulus includes: By optimizing the complex modulus of the rational function approximation jet flow mechanism using a genetic algorithm, the optimal coefficients of the rational function with smaller error compared to the theoretical value in the dominant jet flow frequency band are obtained. The effective fluid modulus is calculated based on the rational function corresponding to the optimal coefficients, and the bulk modulus and shear modulus are obtained by using the effective fluid modulus. The parameters of the complex modulus of the GZM model are optimized by using a genetic algorithm to minimize the error between the complex modulus of the GZM model and the bulk modulus or shear modulus, thereby obtaining the parameters of the complex modulus of the GZM model corresponding to the bulk modulus or shear modulus. Based on the parameters of the complex modulus of the optimized GZM model, the GZM model representations of the bulk modulus and the shear modulus are obtained.

5. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 1, characterized in that, The Zener model is used to characterize the intrinsic viscoelasticity of the solid skeleton, expressed as: = ,in, Indicates the strain relaxation time of the skeleton. For skeleton stress relaxation time, Angular frequency, The imaginary unit, The Zener model factor refers to the intrinsic viscoelasticity of a solid skeleton. ； ； For skeleton quality factor, The center frequency.

6. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 5, characterized in that, The total bulk modulus is the product of the Zener model factor of the intrinsic viscoelasticity of the solid skeleton and the GZM model representation of the bulk modulus; The total shear modulus is the product of the Zener model factor of the intrinsic viscoelasticity of the solid skeleton and the GZM model representation of the shear modulus.

7. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 1, characterized in that, The memory variable is represented as: , It is a time variable. It is a deformation rate function that varies with time. For the first The relaxation strength of each variable, It is the first Time derivatives of each variable =3, For the first The stress relaxation time of each variable.

8. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 7, characterized in that, The first-order differential equation of the memory variable is: The current memory variable is updated by comparing the current strain with the memory variable value of the previous time step. Based on the first-order differential equation, the values ​​of all memory variables are updated at each time step. The updated memory variables are substituted into the time-domain velocity-stress equation system containing all memory variables to calculate the stress and fluid pressure at the current time step. Then, the calculated stress and fluid pressure are used to update the velocity and displacement fields of the solid and fluid to obtain the strain field at the next time step. This process is repeated until the simulation is complete.

Images