Porous rock seismic wave propagation characteristic simulation method based on multi-mechanism coupling

Abstract

The invention belongs to the technical field of seismic exploration, and relates to a porous rock seismic wave propagation characteristic simulation method based on multi-mechanism coupling, which comprises the following steps: simulating the propagation of seismic waves in viscoelastic porous rock by adopting a seismic wave propagation model; wherein the modeling process of the seismic wave propagation model comprises the following steps: coupling a jet flow mechanism with a Biot theory framework; the method comprises the following steps: representing the inherent viscoelasticity of a solid skeleton by adopting a Zener model, and coupling with a jet flow mechanism and a Biot theoretical framework to obtain a total volume modulus and a total shear modulus; calculating a memory variable according to the total volume modulus and the total shear modulus; obtaining a time domain relaxation function containing a memory variable; substituting the time domain relaxation function containing the memory variables into the frequency domain stress-strain relation to obtain a time domain speed-stress equation set containing all the memory variables; and obtaining a seismic wave propagation model by combining the first-order differential equation of the memory variable. According to the invention, wave field characteristics such as reflection and diffraction in a complex heterogeneous medium can be accurately and comprehensively simulated.

Description

Technical Field

[0001] This application belongs to the field of seismic exploration technology, specifically a method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling. Background Technology

[0002] In the field of geophysical exploration, the propagation characteristics of seismic waves in porous and dissipative media are a key foundation for understanding underground structures and resource distribution. When seismic waves propagate in porous media such as fluid-saturated rock, they exhibit significant viscoelastic properties, specifically frequency-dependent energy loss (attenuation) and velocity dispersion. These characteristics are of great significance for many geoscience applications, including environmental monitoring, geothermal exploration, oil and gas extraction, geotechnical engineering, and earthquake research.

[0003] Fluid-saturated rocks are often modeled as a two-phase system of "solid skeleton-pore fluid," in which fluid fills the pore spaces of the solid skeleton. Biot's classic theory was the first to systematically describe the propagation of seismic waves in such a two-phase system, clarifying that the flow of fluid relative to the solid skeleton is one of the core 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 significant across a wide frequency range, and the assumption of uniform porosity in Biot theory is insufficient to characterize the mechanical response of complex pore structures in real rocks. Wave-induced fluid flow, as a key mechanism for viscoelastic attenuation in fluid-saturated rocks, plays a role at multiple scales, with microscale jetting being particularly important. This mechanism originates from fluid exchange between pores with different compressibility or permeability; although compliant pores account for only a small portion of the total pore space, they significantly affect the overall mechanical properties of the rock. For the jetting mechanism, simplified and efficient analytical models have been developed by introducing jet flow aspect ratio parameters to assess the seismic dispersion and attenuation characteristics of isotropic porous rocks. However, fluid-related dissipation mechanisms cannot fully explain energy loss in real rocks. The solid skeleton of rocks is not an ideal elastic body; its intrinsic viscoelastic properties also significantly affect seismic wave attenuation and dispersion, and existing models often fail to fully integrate this key mechanism. Summary of the Invention

[0005] This application provides a method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling. This method solves the problem that it is impossible to capture the coupling effect of three key attenuation mechanisms in a wide frequency range within a unified framework, resulting in incomplete simulation of attenuation peaks and dispersion characteristics, which makes it difficult to meet the needs of complex geological scenarios.

[0006] A method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling, according to an embodiment of this application, 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 was obtained; 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. The total bulk modulus and total shear modulus are obtained. The Zener model is used to characterize the intrinsic viscoelasticity of the solid skeleton. The model is coupled with the GZM model representations of the bulk modulus and the shear modulus 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 to obtain the time-domain relaxation function containing memory variables. 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, the time-domain stress-strain relationship containing memory variables is obtained. Differentiating the time-domain stress-strain relationship with respect to time yields the time-domain velocity-stress equations containing all memory variables. Combining these with the first-order differential equations of the memory variables, the seismic wave propagation model is obtained.

[0007] Furthermore, by coupling the jet flow mechanism with the Biot theoretical framework, the frequency domain stress-strain relationship was 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; 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, This is the effective stress coefficient.

[0008] Furthermore, 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.

[0009] Furthermore, 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, , For effective fluid modulus, The imaginary unit, Angular frequency, For fluid shear viscosity, The bulk modulus of a fluid. This is the complex modulus of the jet flow mechanism.

[0010] Furthermore, the frequency-dependent complex modulus in the frequency domain stress-strain relationship is transformed into a complex modulus GZM model representation, including: 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 using the rational function corresponding to the optimal coefficients, and the bulk modulus and shear modulus are then calculated 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.

[0011] Furthermore, 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.

[0012] Furthermore, 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 expressed in the GZM model form as the Zener model factor of the intrinsic viscoelasticity of the solid skeleton and the shear modulus.

[0013] Furthermore, 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.

[0014] 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.

[0015] 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

[0016] Figure 1 A flowchart illustrating the method provided in the embodiments of this application; 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; 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

[0017] 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.

[0018] 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; A seismic wave propagation model was used to simulate the propagation of seismic waves in viscoelastic porous rocks; The modeling process of the seismic wave propagation model includes: By coupling the jet flow mechanism with the Biot theoretical framework, the frequency domain stress-strain relationship was obtained; 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. 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. 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 to obtain the time-domain relaxation function containing memory variables. 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.

[0019] 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.

[0020] 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.

[0021] 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.

[0022] 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: , , 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.

[0023] 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: ; ; 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; 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: , , , , 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; 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; 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; The frequency-dependent characteristic represents the modulus of the correlation strength between fluid pressure and deformation at different frequencies; 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.

[0024] 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.

[0025] 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. 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. 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: 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.

[0026] 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 ; 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. 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 between the frequency-varying modulus and the jet flow is minimal. The frequency-varying modulus includes the bulk modulus and the shear modulus.

[0027] 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; 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: = = , It is a deformation rate function that varies with time. It is the time derivative of the memory variable.

[0028] 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.

[0029] 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. .

[0030] 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: , , , 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 , 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; 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; 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; 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; 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; 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; 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; It is the first Each Zener element corresponds to a coupling modulus Memory variables; 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; It is the first Effective stress coefficient of each Zener element Memory variables; It is the first Zener element The memory variable of the modulus.

[0031] 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: 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: 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.

[0032] 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: , .

[0033] 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). 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: = = , 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.

[0034] 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.

[0035] 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.

[0036] 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: , , , 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; 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; 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; 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; 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; 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; 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; 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; 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.

[0037] 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.

[0038] 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 which Figure 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. 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 was obtained; 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, wherein 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. 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 to obtain the time-domain relaxation function containing memory variables 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, the time-domain stress-strain relationship containing memory variables is obtained. Differentiating the time-domain stress-strain relationship with respect to time yields the time-domain velocity-stress equations containing all memory variables. Combining these with the first-order differential equations of the memory variables, 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, 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, This is the effective stress coefficient.

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 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.

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, 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.

5. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 4, 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 using the rational function corresponding to the optimal coefficients, and the bulk modulus and shear modulus are then calculated 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.

6. 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.

7. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 6, 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 expressed in the GZM model form as the Zener model factor of the intrinsic viscoelasticity of the solid skeleton and the shear modulus.

8. 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.

9. The method for simulating the propagation characteristics of seismic waves in porous rocks based on multi-mechanism coupling according to claim 8, 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