Fast broadband visco-elasticity forward modeling method for visco-elasticity of vibroseis

Abstract

The invention provides a fast broadband vibroseis viscoelasticity forward modeling method. The fast broadband vibroseis viscoelasticity forward modeling method comprises the steps that 1, a pulse vertical force source function is loaded at a vibroseis excitation point position; 2, calculating a visco-elastic equation Green function of a visco-elastic excitation point location; 3, generating a vibroseis record excited at the position according to the convolution broadband vibroseis function and the viscoelastic equation Green function; and 4, linearly combining the single-shot synthetic seismic records according to an observation system to form a combined-shot vibroseis synthetic seismic record. According to the rapid broadband vibroseis viscoelasticity forward modeling method, the calculation cost of recording broadband vibroseis forward modeling for a long time is effectively reduced, and high-efficiency broadband vibroseis viscoelasticity forward modeling is realized. According to the invention, high-precision and high-efficiency seismic simulation of the vibroseis under complex geological conditions is realized, and an important tool is provided for design and analysis of a vibroseis observation system.

Description

Technical Field

[0001] The invention relates to the technical field of geophysical exploration seismic data acquisition and processing, and in particular to a fast broadband controllable viscoelastic forward modeling method. Background Art

[0002] Viscoacoustic forward modeling technology with controlled seismic source is an important technology in the field of geophysical exploration and is widely used in oil and gas resource exploration, geological disaster warning, and marine resource development.

[0003] In actual seismic exploration, the underground medium has viscosity, and there are attenuation and dispersion effects, which cause the main frequency of seismic waves to decrease, the amplitude to decrease, and the resolution of seismic imaging to decrease, making the underground structure positioning inaccurate and affecting the reliability of seismic interpretation. The viscoelastic wave equation provides a theoretical tool for seismological research. In the viscoelastic wave equation, the viscoelastic properties of the formation need to be considered to calculate the Green's Function, which means that the elastic parameters and viscosity parameters of the formation will affect the calculation results of the Green's Function. Conventional processing methods include finite difference method, finite element method, optical source method, and compressed sensing-based methods. Among them, the finite difference method is widely used in the numerical calculation of Green's function because of its simplicity, stability, and high efficiency. In the field of seismology, the Green's function calculation method of the viscoelastic wave equation and the viscoelastic forward simulation technology of the earthquake source have become an important part of the research.

[0004] The synthesis of seismic records is an important part of seismological research, which can help us understand the propagation of seismic waves and the characteristics of the strata. In practice, the linear time-invariant property of the viscoelastic wave equation can be used to synthesize seismic records by convolution of the source time function and the Green's function. In this process, the source time function describes the change of the intensity of the source over time, while the Green's function describes the propagation process of seismic waves in the strata.

[0005] Traditional narrow-band forward modeling technology only considers the energy transfer on the wave propagation path, but ignores the energy transfer on the side of the wave propagation path. Therefore, it can only simulate seismic records within a narrow frequency band, but cannot simulate wide-band seismic records. It is necessary to perform multiple forward simulations to obtain seismic records within different frequency bands, which will result in a large amount of data and increase the difficulty of data processing and storage. Wide-band forward modeling technology can simulate a wider frequency band and can more comprehensively reflect the characteristics of underground media. Only one forward simulation is needed to obtain seismic records within a wide frequency band, which can greatly reduce the amount of data, make data processing and storage more convenient, and thus improve the detection accuracy of seismic exploration and ensure the reliability of exploration results.

[0006] In recent years, broadband controllable vibroseis viscoelastic forward modeling technology has been widely studied and applied. This technology can simulate seismic waves within a wide frequency range, so it can better simulate the propagation process of seismic waves. At the same time, this technology can also simulate different vibration conditions by controlling the parameters of the seismic source.

[0007] In the Chinese patent application with application number: CN201210408060.8, a numerical simulation method for forward modeling of synchronous scanning wavefield of controllable vibroseis is involved, which belongs to the field of petroleum exploration. The method comprises: (1) inputting simulation parameters; (2) calculating and obtaining the simulation parameters of synchronous scanning of controllable vibroseis according to the simulation parameters input in step (1), and outputting these parameters into a parameter file; the simulation parameters of synchronous scanning of controllable vibroseis include: the position of each source, the start and end scanning time and which synchronous source it belongs to, the continuous recording time length and the scanning signal of each synchronous element; (3) obtaining the simulation record by high-order finite difference wavefield extension calculation; (4) outputting the simulation result. The inventive method can simulate the synchronous scanning of controllable vibroseis in various ways such as sliding scanning, independent synchronous scanning, pseudo-random scanning, distance-pull synchronous scanning technology and V1, and the forward modeling simulation continuous recording of any number of controllable vibroseis according to the needs of actual simulation through the controllable input parameters.

[0008] In the Chinese patent application with application number: CN201510272589.5, a parallel forward modeling method based on a controllable seismic source sliding scanning method is involved, and the parallel forward modeling method based on a controllable seismic source sliding scanning method includes: step 1, input simulation parameters, including conventional observation system parameters and controllable seismic source sliding scanning parameters; step 2, calculate the controllable seismic source sliding scanning simulation parameters according to the input simulation parameters; step 3, perform controllable seismic source sliding scanning forward modeling according to the controllable seismic source sliding scanning simulation parameters; step 4, output the parallel forward modeling results. Compared with conventional serial forward modeling, the numerical simulation results of the parallel forward modeling method based on the controllable seismic source sliding scanning method are more accurate, and the biggest advantage is that the computational efficiency of the present invention is greatly improved compared with the serial forward modeling.

[0009] In the Chinese patent application with application number: CN202111271143.2, a controllable vibroseis forward modeling method, device, storage medium and electronic equipment are involved, including: establishing a velocity model according to the geological characteristics of the target desert area; setting the observation system, source wavelet type, frequency scanning mode, frequency scanning range and signal scanning length according to the velocity model; starting the forward modeling to obtain the scanning signal and the source record; cross-correlating the scanning signal and the source record to obtain a single shot record. This invention can accurately simulate the acquisition of controllable vibroseis data on the desert surface, meet the current research needs of noise suppression of data collected by controllable vibroseis on the desert surface, and provide technical support for the analysis of the formation mechanism of noise in seismic exploration in desert areas and for the targeted suppression of noise.

[0010] The above existing technologies are greatly different from the present invention and fail to solve the technical problem we want to solve. For this reason, we have invented a new fast broadband controllable viscoelastic forward modeling method. Summary of the invention

[0011] The purpose of the present invention is to provide a method for rapidly performing broadband vibroseis viscoelastic forward modeling, thereby providing an important tool for the design and analysis of vibroseis observation systems under complex geological conditions.

[0012] The object of the present invention can be achieved by the following technical measures: a fast broadband controllable viscoelastic forward modeling method, the fast broadband controllable viscoelastic forward modeling method comprising:

[0013] Step 1: Load the pulse vertical force source function at the vibrator excitation point;

[0014] Step 2: Calculate the Green's function of the viscoelastic equation at the vibrator excitation point;

[0015] Step 3: Generate a vibroseis record excited at the location based on the convolution broadband vibroseis function and the Green's function of the viscoelastic equation;

[0016] Step 4: According to the observation system, linearly combine the single-shot synthetic seismic records to form the combined-shot vibroseis synthetic seismic records.

[0017] The purpose of the present invention can also be achieved by the following technical measures:

[0018] In step 1, the pulse vertical force source function is loaded according to the position of the vibrator excitation.

[0019] In step 1, the impulse vertical force source function can be expressed as:

[0020]

[0021] where x 0 =(x 0,y 0 ,z 0 ) is the spatial position of the point force source, g(t,t 0 ,ω) is the source time function, which is the time when it reaches t 0 and frequency ω; (F x ,F y ,F z ) T are the three components of the force vector, F 0 is the amplitude of the force; for a controllable vibrator, only the vertical force is considered, that is, F x =F y = 0, but it is also applicable to source simulations with force vectors in different directions.

[0022] In step 1, the source time function g(t,t 0 ,ω) uses the pulse function δ(tt 0 ); The pulse function is a special function, except that t = t 0 Except for pulses, the function is zero at other times; in addition, the integral of the function is 1; assuming that the function f(t) is a continuous function, the pulse function satisfies the following conditions:

[0023] ∫f(t)δ(tt 0 )dt=f(t 0 ), (2).

[0024] In step 1, at time grid t n =nΔt discretized pulse function, where Δt is also the time step for solving the viscoelastic wave equation; the discretized pulse function sequence is expressed as d n ,n=0,1,2,…; According to the properties of the pulse function, for the polynomial function f(t)=t q ,q=0,1,…,Q,time series d n The following conditions must be met:

[0025]

[0026] According to the Q+1 equations above, we can solve the time series d n ; Solved d n Except for the case where t = t 0 Outside of this neighborhood, there is zero everywhere else.

[0027] In step 2, at the controllable vibrator excitation point, the pulse time function is loaded, the viscoelastic wave equation is solved, and the Green's function of the viscoelastic wave equation is obtained.

[0028] In step 2, a linear viscoelastic model composed of n standard linear solids is considered, and the corresponding control equation is:

[0029]

[0030] Where ρ is the density, u is the displacement vector, λ 0 and μ 0 is the medium parameter at the reference frequency, λ v and μ v is the viscoelastic medium parameter (ν=1,2,…,n), and the spatial operator L is:

[0031]

[0032] Auxiliary variables Solve the following equation to get:

[0033]

[0034] The relaxation frequency ω ν >0,v=1,2,…,n.

[0035] In step 2, the Green's function of the viscoelastic wave equation is calculated, which requires solving 3+3n partial differential equations, where n is the number of standard linear solids; taking into account both accuracy and efficiency, n is generally set to 3.

[0036] In step 3, a single-shot synthetic seismic record can be generated by convolving the time function of the controllable source and the Green's function of the viscoelastic wave equation; since the viscoelastic wave equation is a linear time-invariant system with respect to the source time function, the impulse response, i.e., the Green's function, can fully characterize the characteristics of the forward simulation system of the viscoelastic wave equation.

[0037] In step 3, the vibrator excites a continuous sweep signal, usually a linear sinusoidal sweep signal:

[0038] s(t)=a(t)sin[2πφ(t)+ψ]. (7)

[0039] Where ψ is a constant phase and the amplitude term a(t) is:

[0040]

[0041] Where T is the duration;

[0042] The phase term φ(t) is:

[0043]

[0044] where f 1 ,f 2 They are the starting scanning frequency and the ending scanning frequency respectively.

[0045] The purpose of the present invention can also be achieved through the following technical measures: a fast broadband controllable source viscoelastic forward modeling system, which uses a fast broadband controllable source viscoelastic forward modeling method to perform high-precision and high-efficiency controllable source earthquake simulation under complex geological conditions at the time of oil drilling accidents.

[0046] The fast broadband vibroseis viscoelastic forward modeling method of the present invention numerically calculates the Green's function of the viscoelastic wave equation, utilizes the linear time-invariant characteristics of the viscoelastic wave equation, and synthesizes earthquake records by convolving the source time function. Compared with the existing methods, the present invention has the following advantages: (1) The calculation cost is independent of the length of the source time function, and the calculation efficiency is high; (2) After the Green's function is prepared, it can be applied to different vibroseis observation systems, including different vibroseis spatial distributions, different excitation times, different source time functions, etc., which is very flexible.

[0047] The fast broadband vibroseis viscoelastic forward modeling method of the present invention effectively reduces the computational cost of long-term recording broadband viscoelastic forward simulation, and realizes high-efficiency broadband vibroseis viscoelastic forward simulation. The present invention realizes high-precision and high-efficiency earthquake simulation of vibroseis under complex geological conditions, and provides an important tool for the design and analysis of vibroseis observation systems. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] Figure 1 A flowchart of a specific embodiment of the fast broadband vibroseis viscoelastic forward modeling method of the present invention;

[0049] Figure 2 is a viscoelastic velocity model in a specific embodiment of the present invention;

[0050] Figure 3 is the viscoelastic Green's function calculated in one specific embodiment of the present invention;

[0051] Figure 4 The broadband vibrator time function and its time-frequency spectrum in a specific embodiment of the present invention;

[0052] Figure 5 It is a seismic record synthesized by convolving the viscoelastic Green's function and the vibrator time function in a specific embodiment of the present invention;

[0053] Figure 6 It is the result of cross-correlation between the synthetic record of the vibrator and the time function of the vibrator in a specific embodiment of the present invention;

[0054] Figure 7 The frequency spectrum of the viscoelastic Green's function, the synthetic seismic record and the cross-correlation record in one specific embodiment of the present invention;

[0055] Figure 8It is a synthetic seismic record produced by delayed excitation of a controllable seismic source in a specific embodiment of the present invention. DETAILED DESCRIPTION

[0056] It should be noted that the following detailed descriptions are exemplary and are intended to provide further explanation of the present invention. Unless otherwise specified, all technical and scientific terms used herein have the same meanings as those commonly understood by those skilled in the art to which the present invention belongs.

[0057] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit exemplary embodiments according to the present invention. As used herein, unless the context clearly indicates otherwise, the singular form is also intended to include the plural form. In addition, it should be understood that when the terms "comprising" and / or "including" are used in this specification, it indicates the presence of features, steps, operations and / or combinations thereof.

[0058] In conventional controllable source forward modeling, the source time function is directly loaded to solve the viscoelastic wave equation, and the calculation time depends on the length of the source time function. Since the controllable source is a continuous time signal and may last for a long time, the calculation cost of the existing method will increase linearly with the increase of duration. By calculating the Green's function and utilizing the linear time-invariant characteristics of the forward simulation system, the present invention can effectively avoid the problem of calculation cost. The increase in the length of the controllable source time function will not increase the calculation cost of solving the Green's function, but will only slightly increase the calculation cost of the one-dimensional convolution. Compared with solving the viscoelastic wave equation, the calculation cost of the one-dimensional convolution can be ignored, realizing fast broadband controllable source viscoelastic forward modeling.

[0059] The fast broadband vibroseis viscoelastic forward modeling method of the present invention comprises: numerically calculating the Green's function of the viscoelastic wave equation, utilizing the linear time-invariant characteristics of the viscoelastic wave equation, and synthesizing seismic records by convolving the source time function. The present invention rapidly calculates the synthetic seismic records of the vibroseis under viscoelastic medium conditions, and can be used for the design and analysis of the vibroseis observation system.

[0060] The following are several specific embodiments of the present invention.

[0061] Example 1

[0062] In a specific embodiment 1 of the present invention, Figure 1 As shown, Figure 1 The flowchart of the fast broadband controllable viscoelastic forward modeling method of the present invention is as follows. The fast broadband controllable viscoelastic forward modeling method includes the following steps:

[0063] Step 1: Load the pulse vertical force source function at the vibrator excitation point;

[0064] According to the position of the vibrator excitation, the pulse function is discretized and used as a vertical force source to load the source time function.

[0065] Step 2: Calculate the Green's function of the viscoelastic equation at the vibrator excitation point;

[0066] At the vibrator excitation point, the pulse time function is loaded, the viscoelastic wave equation is solved, and the Green's function of the viscoelastic wave equation is obtained.

[0067] Step 3: Convolve the broadband vibrator function and the Green's function of the viscoelastic equation to generate the vibrator record excited at the location;

[0068] By convolving the time function of the controllable source and the Green's function of the viscoelastic wave equation, a single-shot synthetic seismic record can be generated. Since the viscoelastic wave equation is a linear time-invariant system with respect to the source time function, the impulse response (Green's function) can fully characterize the characteristics of the viscoelastic wave equation forward simulation system.

[0069] Step 4: Combine the vibroseis records at each location to form a combined gun vibroseis record.

[0070] According to the observation system, the single-shot synthetic seismic records can be linearly combined to form the combined-shot controllable source synthetic seismic records.

[0071] Example 2

[0072] In a specific embodiment 2 of the present invention, the fast broadband vibroseis viscoelastic forward modeling method provided by the present invention comprises the following steps:

[0073] 1. Load the pulse vertical force source function at the vibrator excitation point

[0074] According to the position of the vibrator excitation, the pulse vertical force source function is loaded. The pulse vertical force source function can be expressed as:

[0075]

[0076] where x 0 =(x 0 ,y 0 ,z 0 ) is the spatial position of the point force source, g(t,t 0 ,ω) is the source time function, which is the time when it reaches t 0 and frequency ω. (F x ,F y ,F z ) T are the three components of the force vector, F 0 For controllable vibrators, generally only the vertical force is considered, that is, Fx =F y =0, but this method is also applicable to the simulation of earthquake sources with force vectors in different directions.

[0077] The earthquake source time function g(t,t 0 ,ω) uses the pulse function δ(tt 0 ). The impulse function is a special function, except that t = t 0 Except for the pulse, the function is zero at other times; in addition, the integral of the function is 1. Assuming that the function f(t) is a continuous function, the pulse function satisfies the following conditions:

[0078] ∫f(t)δ(tt 0 )dt=f(t 0 ), (2)

[0079] In the time grid t n =nΔt discretizes the impulse function, where Δt is also the time step for solving the viscoelastic wave equation. The discretized impulse function sequence is expressed as d n ,n=0,1,2,… According to the properties of the impulse function, for the polynomial function f(t)=t q ,q=0,1,…,Q,time series d n The following conditions must be met:

[0080]

[0081] According to the Q+1 equations above, we can solve the time series d n The solution to d n Except for the case where t = t 0 Except for the vicinity, all other places are zero. It should be noted that this method does not require t 0 It falls on the time grid point and has good adaptability to the time grid.

[0082] 2. Calculate the Green's function of the viscoelastic equation at the vibrator excitation point

[0083] Consider a linear viscoelastic model composed of n standard linear solids. The corresponding control equation is:

[0084]

[0085] The spatial operator is:

[0086]

[0087] Auxiliary variables Solve the following equation to get:

[0088]

[0089] The relaxation frequency ω v >0,v=1,2,…,n.

[0090] In this step, in order to obtain the Green's function of the viscoelastic wave equation, 3+3n partial differential equations need to be solved, where n is the number of standard linear solids. Taking into account both accuracy and efficiency, n is generally set to 3. Therefore, compared with elastic wave forward modeling (which requires solving 3 partial differential equations), the computational complexity of viscoelastic wave forward modeling is much greater.

[0091] 3. Convolve the broadband vibrator function and the Green's function of the viscoelastic equation to generate the vibrator record excited at that location

[0092] By convolving the time function of the controllable source and the Green's function of the viscoelastic wave equation, a single-shot synthetic seismic record can be generated. Since the viscoelastic wave equation is a linear time-invariant system with respect to the source time function, the impulse response (Green's function) can fully characterize the characteristics of the viscoelastic wave equation forward simulation system.

[0093] In this step, the vibrator excites a continuous scanning signal, usually a linear sinusoidal scanning signal:

[0094] s(t)=a(t)sin[2πφ(t)+ψ]. (7)

[0095] The amplitude term is:

[0096]

[0097] Where T is the duration.

[0098] The phase term is:

[0099]

[0100] where f 1 ,f 2 They are the starting scanning frequency and the ending scanning frequency respectively.

[0101] 4. Combine the vibroseis records at each location to form a combined gun vibroseis record

[0102] By combining single-shot synthetic seismic records, a combined-shot controlled-source synthetic seismic record can be formed. In conventional controlled-source forward modeling, the source time function is directly loaded to solve the viscoelastic wave equation, and the calculation time depends on the length of the source time function. Since the controlled-source is a continuous time signal that may last for a long time, the computational cost of the existing method increases linearly with the increase in duration. By calculating the Green's function and utilizing the linear time-invariant characteristics of the forward modeling system, the present invention can effectively avoid the computational cost problem. The increase in the length of the controlled-source time function will not increase the computational cost of solving the Green's function, but will only slightly increase the computational cost of the one-dimensional convolution. Compared with solving the viscoelastic wave equation, the computational cost of the one-dimensional convolution can be ignored, realizing fast broadband controlled-source viscoelastic forward modeling.

[0103] Example 3

[0104] In a specific embodiment 3 of the present invention, the fast broadband vibroseis viscoelastic forward modeling method of the present invention is described in detail:

[0105] Figure 2 The viscoelastic P-wave velocity model used in this test is designed based on real geological conditions and contains various complex geological features, such as thin interlayers, pinch-outs, and strong reflection interfaces. The S-wave velocity, density, P-wave quality factor, and S-wave quality factor are obtained by empirical formulas. Figure 3 The viscoelastic Green's function obtained by loading a pulse vertical force source at the middle position of the model surface (x = 4000m). It can be seen from the figure that the viscoelastic seismic wave field is very complex, including strong surface waves and surface wave scattering energy. In addition, affected by the viscoelastic medium itself, the reflected wave energy decays rapidly. This is bound to bring challenges to the imaging of deep reflective layers.

[0106] Figure 4 This is a common vibrator sweep signal and its amplitude spectrum. The starting and ending frequencies are 1 Hz and 60 Hz respectively. Figure 5 is the synthetic seismic record obtained by convolving the vibrator time function and the Green function. Since the vibrator time function is a duration signal, the seismic record after convolution also has an obvious tailing effect. Generally, before processing the vibrator seismic record, it is necessary to cross-correlate the vibrator record and the source time function. The result after cross-correlation is as follows Figure 6 As shown in the figure, it can be seen that the tailing effect of the event axis is suppressed. Figure 7 They are respectively the amplitude spectra of Green's function, vibroseis synthetic seismic record, and cross-correlation result. As can be seen from the figure, the convolution and cross-correlation operations basically do not destroy the characteristics of the amplitude spectrum, thereby verifying the correctness and reliability of the present invention. The present invention is further applied to different vibroseis observation methods. Figure 8Two different controllable vibroseis observation methods are shown. If the conventional controllable vibroseis viscoelastic forward modeling method is used, three 3.5-second viscoelastic forward modeling calculations are required, while the method introduced in the present invention only requires two 3.5-second viscoelastic forward modeling calculations, which greatly reduces the calculation cost.

[0107] Finally, it should be noted that the above is only a preferred embodiment of the present invention and is not intended to limit the present invention. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art can still modify the technical solutions recorded in the aforementioned embodiments or replace some of the technical features therein with equivalents. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

[0108] Except for the technical features described in the specification, all other technical features are known technologies to those skilled in the art.

Claims

1. A fast broadband vibroseis viscoelastic forward modeling method, characterized in that: The fast broadband vibroseis viscoelastic forward modeling method includes: Step 1: Load the pulse vertical force source function at the vibrator excitation point; Step 2: Calculate the Green's function of the viscoelastic equation at the vibrator excitation point; Step 3: Generate a vibroseis record excited at the location based on the convolution broadband vibroseis function and the Green's function of the viscoelastic equation; Step 4: According to the observation system, linearly combine the single-shot synthetic seismic records to form the combined-shot vibroseis synthetic seismic records.

2. The fast broadband vibroseis viscoelastic forward modeling method according to claim 1, characterized in that: In step 1, the pulse vertical force source function is loaded according to the position of the vibrator excitation.

3. The fast broadband vibroseis viscoelastic forward modeling method according to claim 2, characterized in that: In step 1, the impulse vertical force source function can be expressed as: Where x0 = (x0, y0, z0) is the spatial position of the point force source, g(t, t0, ω) is the source time function, which is a function of the arrival time t0 and the frequency ω; (F x ,F y ,F z ) T are the three components of the force vector, F0 is the amplitude of the force; for the controllable vibrator, only the vertical force is considered, that is, F x =F y = 0, but it is also applicable to source simulations with force vectors in different directions.

4. The fast broadband vibroseis viscoelastic forward modeling method according to claim 3, characterized in that: In step 1, the earthquake source time function g(t, t0, ω) adopts the pulse function δ(t-t0); the pulse function is a special function, which is zero at all times except t=t0; in addition, the integral of the function is 1; assuming that the function f(t) is a continuous function, the pulse function satisfies the following conditions: ∫f(t)δ(t-t0)dt=f(t0), (2).

5. The fast broadband vibroseis viscoelastic forward modeling method according to claim 4, characterized in that: In step 1, at time grid t n =nΔt discretized pulse function, where Δt is also the time step for solving the viscoelastic wave equation; the discretized pulse function sequence is expressed as d n ,n=0,1,2,…; According to the properties of the pulse function, for the polynomial function f(t)=t q ,q=0,1,…,Q,time series d n The following conditions must be met: According to the Q+1 equations above, we can solve the time series d n ; Solved d n Except near t=t0, it is zero everywhere else.

6. The fast broadband vibroseis viscoelastic forward modeling method according to claim 1, characterized in that: In step 2, at the controllable vibrator excitation point, the pulse time function is loaded, the viscoelastic wave equation is solved, and the Green's function of the viscoelastic wave equation is obtained.

7. The fast broadband vibroseis viscoelastic forward modeling method according to claim 6, characterized in that: In step 2, a linear viscoelastic model composed of n standard linear solids is considered, and the corresponding control equation is: Where ρ is the density, u is the displacement vector, λ0 and μ0 are the medium parameters at the reference frequency, λ ν and μ ν is the viscoelastic medium parameter (v=1,2,…,n), and the spatial operator L is: Auxiliary variables Solve the following equation to get: The relaxation frequency ω v >0,v=1,2,…,n.

8. The fast broadband vibroseis viscoelastic forward modeling method according to claim 7, characterized in that: In step 2, the Green's function of the viscoelastic wave equation is calculated, which requires solving 3+3n partial differential equations, where n is the number of standard linear solids; taking into account both accuracy and efficiency, n is generally set to 3.

9. The fast broadband vibroseis viscoelastic forward modeling method according to claim 1, characterized in that: In step 3, a single-shot synthetic seismic record can be generated by convolving the time function of the controllable source and the Green's function of the viscoelastic wave equation; since the viscoelastic wave equation is a linear time-invariant system with respect to the source time function, the impulse response, i.e., the Green's function, can fully characterize the characteristics of the forward simulation system of the viscoelastic wave equation.

10. The fast broadband vibroseis viscoelastic forward modeling method according to claim 9, characterized in that: In step 3, the vibrator excites a continuous sweep signal, usually a linear sinusoidal sweep signal: s(t)=a(t)sin[2πφ(t)+ψ]. (7) Where ψ is a constant phase and the amplitude term a(t) is: Where T is the duration; The phase term φ(t) is: Where f1 and f2 are the starting scanning frequency and the ending scanning frequency respectively.

11. Fast broadband controllable viscoelastic forward modeling system, characterized by: The fast broadband vibroseis viscoelastic forward modeling system adopts the fast broadband vibroseis viscoelastic forward modeling method described in any one of claims 1 to 10 to perform high-precision and high-efficiency vibroseis seismic simulation under complex geological conditions for petroleum seismic exploration.

Images