New elastic wave least-squares reverse-time migration method, system, device, and medium

Abstract

The present application relates to a kind of new elastic wave least square inverse time migration method, system, device and medium, comprising: obtaining longitudinal wave velocity model, transverse wave velocity model and density model and observation system OBS 4C data;Wave field is based on longitudinal wave velocity model, transverse wave velocity model and density model and is forward extended;Residual reverse extension is carried out based on longitudinal wave velocity model, transverse wave velocity model and density model and 4C data;Forward extension and residual reverse extension are correlated to obtain initial multi-parameter imaging;Initial multi-parameter imaging is based on ACE equation Born forward simulation and obtains OBS 4C data, the OBS 4C data obtained by simulation is with observation system OBS 4C data in residual reverse extension, and forward extension and residual reverse extension are correlated to calculate, after multiple iterations, LSM imaging result is obtained, and then seismic imaging result is obtained.

Description

Technical Field

[0001] This invention relates to a novel elastic wave least squares reverse time migration method, system, equipment, and medium based on a modified fluctuating acoustic-elastic coupling equation, and pertains to the field of petroleum geophysical exploration technology. Background Technology

[0002] With the rapid development of data acquisition methods and processing technologies in modern offshore oil and gas seismic exploration, offshore oil and gas seismic exploration is evolving from shallow water to deep water, from sea surface to seabed, and from single-component to multi-component exploration. In recent years, deep-water oil and gas seismic exploration has become a key focus and challenge in exploration geophysics both domestically and internationally.

[0003] Subsea regions are elastic, and single-component acquisition by towed cables cannot meet the needs of the oil and gas industry. Submarine seismic exploration deploys four-component (4C) sensors, including three-component (3C) geophones and hydrophones on the seabed interface, and quantitatively estimates subsurface elastic parameters by processing the observed seabed 4C data. Wide azimuth, multi-coverage, high signal-to-noise ratio, and multi-component seabed observation have become important geophysical methods for offshore oil and gas exploration, development, and estimation. The greatest advantage of seabed seismic acquisition is that the 4C seismic data observed above the seabed carries effective subsurface elastic information. Therefore, this multi-component data can be used for elastic wave imaging.

[0004] Traditional acoustic wave equations in marine environments lack information about shear waves, while single elastic wave equations are unstable and prone to numerical dispersion. In addition, rugged seafloor interfaces must be overcome, as their dramatic undulations severely interfere with imaging results. These complex boundary structures generate strong scattering, leading to complex wave fields and waveform distortions, which significantly complicates subsequent seismic processing. Summary of the Invention

[0005] This invention aims to at least solve one of the technical problems existing in the prior art. Therefore, in response to the above-mentioned problems, the purpose of this invention is to provide a novel least-squares reverse-time migration method for elastic waves based on a modified undulating acoustic-elastic coupling equation. This method can solve the problems of low computational efficiency in acoustic-elastic coupling media, artificial corner scattering caused by boundaries, and the lack of shear wave information and instability in conventional single wave equations, thus obtaining accurate imaging results.

[0006] To achieve the aforementioned objectives, the technical solution of this invention is: a novel least-squares reverse-time migration method for elastic waves, characterized by comprising:

[0007] Acquire P-wave velocity model, S-wave velocity model, and density model, as well as data from the OBS 4C observation system;

[0008] Forward wavefield extension based on P-wave velocity model, S-wave velocity model and density model;

[0009] Residual inverse extrapolation is performed based on the P-wave velocity model, S-wave velocity model, density model, and 4C data;

[0010] The initial multi-parameter imaging is obtained by cross-correlation calculation of forward continuation and residual backward continuation;

[0011] The initial multi-parameter imaging was based on the Born forward modeling of the ACE equation to obtain OBS 4C data. The simulated OBS 4C data and the observation system OBS 4C data were subjected to residual backward continuation. The forward continuation and the residual backward continuation were cross-correlated and calculated. After multiple iterations, the LSM imaging results were obtained, and then the seismic imaging results were obtained.

[0012] Furthermore, the forward wavefield extension based on the P-wave velocity model, S-wave velocity model, and density model is performed using the following formula:

[0013]

[0014]

[0015]

[0016] Where (x,z) are the Cartesian coordinates of the grid points before the coordinate transformation, and (ξ,η) are the coordinates of the grid points after the coordinate transformation; and Here, ρ is the directional partial derivative, ρ0 is the density, t is the time, f is the source term, and τ is the source term. xx0 and τ zz0 τ represents the normal stress component. xz0 Represents the shear stress component, v x0 v z0 λ0 and μ0 represent the velocity components in the horizontal x-direction and vertical z-direction, respectively. P0 is the acoustic stress, and λ0 and μ0 are Lamé constants.

[0017] Furthermore, the calculation formula for residual inverse extrapolation based on the P-wave velocity mode, S-wave velocity model, density model, and 4C data is as follows:

[0018]

[0019] Wherein, the superscript R represents the wave field value at the receiver location, x r Indicates the location of the detector point, Δd x (x r ,t) and Δd z (x r ,t) represent the composite data d in the horizontal and vertical directions, respectively. calx With observation data d obsxThe difference between them, (x,z) are the Cartesian coordinates of the grid points before the coordinate transformation, (ξ,η) are the body-fitted coordinates of the grid points after the coordinate transformation, t is time, and τ is the time. xx and τ zz τ represents the normal stress component; xz Represents the shear stress component, v x v z λ and μ represent the velocity components in the horizontal x-direction and vertical z-direction, respectively, and are Lamé constants.

[0020] Furthermore, the initial multi-parameter imaging, based on the ACE equation Born forward modeling, yields OBS 4C data using the following formula:

[0021]

[0022]

[0023]

[0024] Where δ represents the perturbation, (x,z) are the Cartesian coordinates of the grid points before coordinate transformation, (ξ,η) are the body-fitted coordinates of the grid points after coordinate transformation, P is the acoustic stress, ρ0 is the density, t is the time, τxx and τzz represent the normal stress components; τxz represents the shear stress components, vx0 and vz0 represent the velocity components in the horizontal x-direction and vertical z-direction, respectively, λ0 and μ0 are Lamé constants, and m(λ)=δλ and m(μ)=δμ are the perturbations of λ and μ, respectively.

[0025] Furthermore, the forward continuation and the residual backward continuation are cross-correlated and iterated multiple times to obtain the LSM imaging results, which are then output as seismic imaging results, including:

[0026] Construct an objective function and determine whether the objective function meets the set requirements. If it does, output the final seismic imaging result. Otherwise, apply the diagonal Hessian operator to update the LSM imaging result. Calculate the step size using the conjugate gradient method and perform multiple iterations until the set requirements are met, then output the final seismic imaging result.

[0027] Furthermore, the gradient formula for updating the LSM imaging results is as follows:

[0028]

[0029]

[0030] Where g represents the gradient direction, λ and μ are Lamé constants, (x0, z0) are the Cartesian coordinates of the grid points before coordinate transformation, (ξ, η) are the body-fitted coordinates of the grid points after coordinate transformation, and v x vz These represent the velocity components in the horizontal x-direction and vertical z-direction, respectively. The superscript R indicates the wave field value at the receiver location, and τ xx and τ zz τ represents the normal stress component. xz Let represent the shear stress component, t be time, and E be the objective function.

[0031] Secondly, this invention provides a novel elastic wave least-squares reverse-time migration system, the system comprising:

[0032] The first processing unit is configured to acquire P-wave velocity model, S-wave velocity model and density model, as well as OBS 4C observation system data;

[0033] The second processing unit is configured to perform forward wavefield extension based on the longitudinal wave velocity model, the transverse wave velocity model, and the density model.

[0034] The third processing unit is configured to perform residual inverse extrapolation based on the P-wave velocity model, S-wave velocity model, density model, and 4C data.

[0035] The fourth processing unit is configured to perform cross-correlation calculations on the forward continuation and the residual backward continuation to obtain the initial multi-parameter imaging;

[0036] The fifth processing unit is configured to obtain OBS4C data from the initial multi-parameter imaging based on the Born forward modeling of the ACE equation, perform reverse data residual extrapolation on the simulated OBS4C data and the OBS4C data from the observation system, and perform cross-correlation calculation on the forward extrapolation and the reverse data residual extrapolation. After multiple iterations, the LSM imaging results are obtained, and then the seismic imaging results are obtained.

[0037] Thirdly, the present invention provides a computer-readable storage medium for storing one or more programs, characterized in that the one or more programs include instructions that, when executed by a computing device, cause the computing device to perform any of the methods.

[0038] Fourthly, the present invention provides an electronic device comprising: one or more processors, a memory, and one or more programs, wherein the one or more programs are stored in the memory and configured to be executed by the one or more processors, and the one or more programs include instructions for performing any of the methods.

[0039] Because the present invention adopts the above technical solution, it has the following characteristics:

[0040] 1. This invention greatly improves the accuracy of acoustic-elastic coupling medium imaging results and accelerates the convergence speed of the target functional, providing an effective guarantee for reservoir prediction and effectively improving the resolution of seismic data.

[0041] 2. This invention can solve the problems of low computational efficiency of acoustic-elastic coupling media, artificial corner scattering caused by boundaries, and lack of transverse wave information and instability of conventional single wave equations, thus obtaining accurate imaging results.

[0042] 3. This invention utilizes the modified AEC equation to perform a single wavefield simulation prediction and migration of OBS 4C data. This method can generate multi-component images and quantitatively estimate subsurface impedance reflectivity with a computational load similar to the conventional 3C ELSRTM method. Numerical experiments on seabed models and actual models show that this method can improve amplitude fidelity, increase spatial resolution, and suppress migration noise.

[0043] In summary, this invention can be widely applied to offshore oil and gas seismic exploration. Attached Figure Description

[0044] Various other advantages and benefits will become apparent to those skilled in the art upon reading the following detailed description of preferred embodiments. The accompanying drawings are for illustrative purposes only and are not intended to limit the invention. Throughout the drawings, the same reference numerals denote the same parts. In the drawings:

[0045] Figure 1 This is a flowchart of a novel elastic wave least-squares reverse time migration method based on a modified undulating acoustic-elastic coupling equation, according to an embodiment of the present invention.

[0046] Figure 2 This is a simplified undulating seabed model according to an embodiment of the present invention;

[0047] Figure 3 The image shows a 0.6s wavefield snapshot obtained from a simple undulating seabed model in an embodiment of the present invention. (a) and (d) are the acoustic-elastic coupling equations in curved coordinates; (b) and (e) are the elastic wave equations in curved coordinates; (c) and (f) are the acoustic-elastic coupling equations in rectangular coordinates; (a), (b), and (c) are the horizontal components; (d), (e), and (f) are the longitudinal components.

[0048] Figure 4 This is an example of an underwater mountain-valley model in an embodiment of the present invention;

[0049] Figure 5 The following is a body-fitted mesh subdivision of the underwater mountain-valley model in an embodiment of the present invention: (a) is the whole, and (b) is a local magnification.

[0050] Figure 6The images shown are the results of 10 iterations of least-squares reverse-time migration imaging of elastic waves in this embodiment of the invention. (a) is the overall image, and (b) is a magnified local image.

[0051] Figure 7 This is a structural diagram of an electronic device according to an embodiment of the present invention. Detailed Implementation

[0052] It should be understood that the terminology used herein is for the purpose of describing particular exemplary embodiments only and is not intended to be limiting. Unless the context clearly indicates otherwise, the singular forms “a,” “an,” and “described” as used herein may also include the plural forms. The terms “comprising,” “including,” “containing,” and “having” are inclusive and therefore indicate the presence of the stated features, steps, operations, elements, and / or components, but do not exclude the presence or addition of one or more other features, steps, operations, elements, components, and / or combinations thereof. The method steps, processes, and operations described herein are not construed as requiring them to be performed in a particular order described or illustrated unless the order of performance is explicitly indicated. It should also be understood that additional or alternative steps may be used.

[0053] Although terms such as first, second, third, etc., may be used in this document to describe multiple elements, components, regions, layers, and / or segments, these elements, components, regions, layers, and / or segments should not be limited by these terms. These terms may be used only to distinguish one element, component, region, layer, or segment from another. Unless the context clearly indicates otherwise, terms such as "first," "second," and other numerical terms used herein do not imply order or sequence. Therefore, the first element, component, region, layer, or segment discussed below may be referred to as the second element, component, region, layer, or segment without departing from the teachings of the exemplary embodiments.

[0054] For ease of description, spatial relative terms may be used in the text to describe the relationship of one element or feature relative to another element or feature as shown in the figure. These relative terms include, for example, "inside," "outside," "middle," "outer," "below," "above," etc. Such spatial relative terms are intended to include different orientations of the device in use or operation, other than those depicted in the figure.

[0055] Due to the low computational efficiency of existing acoustic-elastic coupling media, artificial corner scattering caused by boundaries, and the lack of shear wave information and instability of conventional single wave equations, this invention provides a novel elastic wave least-squares reverse time migration (ELSRTM) method, system, equipment, and medium based on a modified fluctuating acoustic-elastic coupling (AEC) equation. The method includes: prediction or migration of OBS4C data can be achieved using a single wavefield simulation; combining perturbation imaging conditions to invert P-wave and S-wave velocities and densities, and generating amplitude-preserving imaging results; and implementing a preconditional conjugate gradient algorithm based on a multi-parameter diagonal Hessian matrix to obtain seismic imaging results. This invention significantly improves the accuracy of acoustic-elastic coupling medium imaging results and accelerates the convergence speed of the target functional, providing effective assurance for reservoir prediction and improving the resolution of seismic data.

[0056] Exemplary embodiments of the invention will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the invention are shown in the drawings, it should be understood that the invention can be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided to enable a more thorough understanding of the invention and to fully convey the scope of the invention to those skilled in the art.

[0057] This invention introduces the idea of ​​mapping irregular interfaces to regular interfaces in curved coordinates (in order to adapt to undulating interfaces, all grids are changed to curved grids instead of rectangular grids), derives the modified undulating acoustic-elastic coupling equation in the curved coordinate system, and further extends the equation to least squares reverse time migration, thereby obtaining high-resolution seismic imaging results.

[0058] Based on classical migration theory, the linear Born model can be used to predict seismic data d. obs :

[0059] d obs =Lm (1)

[0060] Where L is the Born operator and m is the reflection coefficient model.

[0061] By separating the parameters and wavefield into background and perturbation components, the elastic Born modeling equation for OBS 4C data can be derived. Parameter perturbations can generate wavefield perturbations, which are as follows:

[0062] δV(λ,μ)=V(λ,μ)-V0(λ,μ) (2)

[0063] δU=U-U0 (3)

[0064] Where U0 and V0 represent the background wave field and background parameter, respectively, δU and δV represent the perturbation wave field and perturbation parameter, respectively, and λ and μ are Lamé constants.

[0065] First, the acoustic-elastic coupling equations in curved coordinates are given for forward wavefield continuation:

[0066]

[0067]

[0068]

[0069] Where (x,z) are the Cartesian coordinates of the grid points before the coordinate transformation, and (ξ,η) are the coordinates of the points in the curved coordinate system after the coordinate transformation. and Here, ρ is the directional partial derivative, ρ0 is the density, t is the time, f is the source term, and τ is the source term. xx0 and τ zz0 τ represents the normal stress component. xz0 Represents the shear stress component, v x0 v z0 λ0 and μ0 represent the velocity components in the horizontal x-direction and vertical z-direction, respectively. P0 is the acoustic stress, and λ0 and μ0 are Lamé constants.

[0070] Dividing the wavefield into background and perturbation components allows us to derive the elastic Born modeling equation for OBS 4C data, thus obtaining the OBS 4C data:

[0071]

[0072]

[0073]

[0074] Where δ represents the perturbation, (x,z) are the Cartesian coordinates of the grid points before coordinate transformation, (ξ,η) are the body-fitted coordinates of the grid points after coordinate transformation, P is the acoustic stress, ρ0 is the density, t is the time, τxx and τzz represent the normal stress components; τxz represents the shear stress components, vx0 and vz0 represent the velocity components in the horizontal x-direction and vertical z-direction, respectively, λ0 and μ0 are Lamé constants, and m(λ)=δλ and m(μ)=δμ are the perturbations of λ and μ, respectively.

[0075] Specifically, a preconditional conjugate gradient algorithm is implemented based on a multi-parameter diagonal Hessian matrix, and a least-squares offset operator in curved coordinates is derived for calculating multi-parameter imaging.

[0076]

[0077] Wherein, the superscript R represents the wave field value at the receiver location, x rIndicates the location of the detector point, Δd x (x r ,t) and Δd z (x r ,t) represent the composite data d in the horizontal and vertical directions, respectively. calx With observation data d obsx The difference between them, (x,z) are the Cartesian coordinates of the grid points before the coordinate transformation, (ξ,η) are the body-fitted coordinates of the grid points after the coordinate transformation, t is time, and τ is the time. xx and τ zz τ represents the normal stress component; xz Represents the shear stress component, v x v z λ and μ represent the velocity components in the horizontal x-direction and vertical z-direction, respectively, and are Lamé constants.

[0078] It can be obtained from the following formula:

[0079]

[0080] Where, x r Indicates the location of the detector point, d calx d represents the composite data in the x-direction. calz d represents the composite data in the z-direction. obsx d represents the observation data in the x-direction. obsz Δd represents the observation data in the z-direction. x Δd represents the data residual in the x-direction. z This represents the data residual in the z-direction, where t is time and X is... r Indicates the location of the detector point.

[0081] To utilize both the observed pressure and displacement components simultaneously, the objective function E(m) from the OBS 4C data was employed. The perturbation δE(m) of the objective function can be expressed as:

[0082]

[0083] Where B represents the operator that restricts the shot record from the entire model space to the detector position, m represents the reflection coefficient, L represents the linear born operator, U represents the wave field value, δ represents the perturbation, and E represents the objective function. From equation (8), we can obtain:

[0084] δU=LF' (9)

[0085] Where L represents the linear born operator, U represents the wave field value, δ represents the perturbation, and F' is the source matrix that generates the perturbation wave field, obtained by the following formula:

[0086]

[0087] Where δ represents the perturbation, λ and μ are Lamé constants, (x0, z0) are the Cartesian coordinates of the grid points before coordinate transformation, (ξ, η) are the body-fitted coordinates of the grid points after coordinate transformation, and F' is the source matrix that generates the perturbation wave field. Then equation (9) can be rewritten as:

[0088] δE(m)= <B(LF'),(BU-d obs )>= <F',L R B R (BU-d obs (11)

[0089] Where F' is the source matrix that generates the perturbation wave field, R represents the conjugate operator, and B R This indicates that the data at the detector point is spread throughout the entire model space, L R Let d represent the inverse time delay topology operator of the wave field, m be the reflection coefficient model, E be the objective function, U be the wave field value, and d be the wave field value. obs Representing the observed data, δE(λ,μ) is:

[0090]

[0091] Where δ represents the perturbation, λ and μ are Lamé constants, (x0, z0) are the Cartesian coordinates of the grid points before the coordinate transformation, (ξ, η) are the body-fitted coordinates of the grid points after the coordinate transformation, and τ xx and τ zz τ represents the normal stress component. xz Represents the shear stress component, v x v z These represent the velocity components in the horizontal x-direction and the vertical z-direction, respectively. The superscript R indicates the wave field value at the receiver location, and E represents the objective function.

[0092] The gradient formula is used to update LSM imaging results. The gradient formula is:

[0093]

[0094]

[0095] Where g represents the gradient direction, λ and μ are Lamé constants, (x0, z0) are the Cartesian coordinates of the grid points before coordinate transformation, (ξ, η) are the body-fitted coordinates of the grid points after coordinate transformation, and v x v z These represent the velocity components in the horizontal x-direction and vertical z-direction, respectively. The superscript R indicates the wave field value at the receiver location, and τ xx and τ zz τ represents the normal stress component. xzLet represent the shear stress component, t be time, and E be the objective function.

[0096] Example 1: Based on the above theoretical derivation, such as Figure 1 As shown, the novel elastic wave least-squares reverse-time migration method based on the modified fluctuating sound-elastic coupling equation provided in this embodiment includes:

[0097] S1, Input longitudinal wave velocity model V P Shear wave velocity model V S And density model and input observation system OBS4C data.

[0098] S2, Based on the longitudinal wave velocity model V P Shear wave velocity model V S The density model uses formula (4) for positive wave field extension, and is based on the longitudinal wave velocity model V. P Shear wave velocity model V S The density model and 4C data are used to perform residual reverse continuation using formula (6), and the forward and reverse continuations are cross-correlated to obtain multi-parameter imaging.

[0099] S3. Perform ACE equation Born forward modeling on multi-parameter imaging and simulate using formula (5) to obtain OBS 4C data;

[0100] S4. The simulated OBS 4C data and the input observation system OBS 4C data are used to perform reverse continuation of the 4C data residuals using formula (6). The cross-correlation of the forward and reverse continuation 4C data residuals is then performed to calculate the LSM imaging results. The high-resolution seismic imaging results are output to reflect the geological structure. Subsequently, by superimposing the imaging results obtained from each shot, the solid region reflection coefficient model is obtained. This model can accurately reflect the geological structure.

[0101] In a preferred embodiment, the specific process of obtaining the LSM imaging result by cross-correlation iterative calculation of the forward and reverse continuation 4C data residuals in step S4 is as follows: construct the objective function E and determine whether the objective function meets the requirements. If it does, output the final imaging result. Otherwise, use formula (13) to apply the diagonal Hessian operator to update the LSM imaging result. Use the conjugate gradient method to obtain the step size and iterate to meet the set requirements, and then update to obtain the LSM imaging result.

[0102] To more intuitively observe the propagation of seismic waves in an acoustic-elastic coupling medium, such as Figure 2 As shown, the following uses a layered model with a horizontal seabed interface as a specific example to numerically simulate the acoustic-elastic coupling equation, the acoustic wave equation, and the elastic wave equation using the second-order time and fourth-order spatial finite difference method, and analyzes the simulation results.

[0103] The model size is 201×151, the spatial grid spacing is Δx=Δz=4m, the time step is Δt=0.3ms, the Ricker wavelet with a dominant frequency of 25Hz is used, the source is placed at (400m, 4m), the acoustic-bullet coupling interface depth is z=160m, and the OBC is located at the seabed horizontal interface.

[0104] from Figure 3 As can be seen in (a), when the sound wave reaches the seabed reflecting interface, it generates ② reflected PP wave and ③ transmitted PP wave. When the transmitted PP wave continues to propagate to the second reflecting interface in the solid medium region, it generates both reflected PP wave and transmitted PP wave again, but no transverse wave information is generated. Figure 3 In (b) and 3(c), seismic waves are excited at sea level and propagate to the seabed interface as sound waves (① direct P-waves). Part of the propagation generates ② reflected P-waves, while the other part is converted into elastic waves, i.e., ③ transmitted P-waves and ④ P-wave converted to S-waves (PS-waves). Since S-waves are slower than P-waves, at T = 330 ms, the P-waves reach the reflection interface in the solid medium first, undergoing emission and transmission, generating ② reflected P-waves, ⑤ reflected PS-waves, ③ transmitted P-waves, and ④ transmitted PS-waves. From the pure elastic wave field snapshots 3(d) and 3(e), it can be seen that when the seabed interface is horizontal, the wave field snapshots of the elastic wave equation and the acoustic-elastic coupling equation are the same. However, in the elastic wave forward modeling, the S-wave velocity in the seawater portion is set to 0, requiring the Poisson's ratio of both P-waves and S-waves to be equal to 0.5 to ensure the stability of the numerical method. Furthermore, setting the S-wave velocity to 0 limits the continuity of velocity in the model. Therefore, this method is not easily implemented when the acoustic-elastic coupling medium is complex. In addition, the computation times for the acoustic wave equation, the acoustic-elastic coupling equation, and the single elastic wave equation were 105.641 s, 133.844 s, and 183.125 s, respectively. This indicates that although the acoustic wave equation requires the shortest computation time, it neglects shear wave information, thus failing to obtain high-resolution imaging results of the solid region below the seabed interface. The elastic wave equation takes even longer to compute than the acoustic-elastic coupling equation, and this difference widens as numerical simulations are performed on larger models or based on real-world data; that is, pure elastic waves require significantly more computation time. Therefore, the acoustic-elastic coupling equation is the most suitable for simulating marine environments.

[0105] In addition, this embodiment tests a seafloor peak-valley model. This model features an undulating seafloor interface; above the interface is a fluid medium, and below it is a solid medium. Figure 4 As can be seen, the model size is 301m × 201m, the grid spacing is Δx = Δz = 8m, and the undulating seabed interface is represented by a function. The seismic source used a Ricker wavelet with a dominant frequency of 30Hz, excited at z = 8m above sea level. A total of 51 shots were used, each with 300 receiver channels, with a shot spacing of 48m and a time step of Δt = 0.5ms. The P-wave velocities were 1500, 2200 m / s, 3100 m / s, and 4000 m / s, respectively; the S-wave velocities were 0, 1270 m / s, 1790 m / s, and 2310 m / s, respectively; and the densities were ρ = 1000 kg / m³. 3 ρ = 2000 kg / m 3 ρ = 2400 kg / m 3 ρ = 2700 kg / m 3 First, the underwater mountain-valley model is discretized, and then a body-fitted mesh is generated by solving the Poisson equation, such as... Figure 5 As shown. To demonstrate the orthogonality of the body-fitted mesh at the undulating seabed boundary, the mountain peaks in the mesh are magnified. Figure 5 As can be seen in (b), its orthogonality at the boundary is very good, the generated mesh is compatible with the undulating acoustic-elastic coupling interface, and the mesh continuity is also very good, which lays the foundation for the stability of subsequent numerical simulations. Figure 6 The results show the PP and SS components of the least squares reverse time migration of the acoustic-elastic coupling equation in the Cartesian coordinate system after inverse transformation of the results in the curved coordinate system. Although the imaging results have low-frequency noise, they basically restore the underground medium structure and do not produce diffraction. Moreover, the deep phase axis is continuous, which shows the correctness of the method in this embodiment.

[0106] Example 2: Following the novel elastic wave least squares reverse time migration method provided in Example 1, this example provides a novel elastic wave least squares reverse time migration system. The system provided in this example can implement the novel elastic wave least squares reverse time migration method of Example 1. This system can be implemented through software, hardware, or a combination of both. For ease of description, this example is described by dividing the system into various functional units. Of course, in implementation, the functions of each unit can be implemented in one or more software and / or hardware components. For example, the system can include integrated or separate functional modules or units to execute the corresponding steps in the methods of Example 1. Since the system in this example is basically similar to the method example, the description process of this example is relatively simple. For relevant details, please refer to the description in Example 1. The example of the novel elastic wave least squares reverse time migration system provided by this invention is merely illustrative.

[0107] This embodiment provides a novel elastic wave least-squares reverse time migration system, which includes:

[0108] The first processing unit is configured to acquire P-wave velocity model, S-wave velocity model and density model, as well as OBS 4C observation system data;

[0109] The second processing unit is configured to perform forward wavefield extension based on the longitudinal wave velocity model, the transverse wave velocity model, and the density model.

[0110] The third processing unit is configured to perform residual inverse extrapolation based on the P-wave velocity model, S-wave velocity model, density model, and 4C data.

[0111] The fourth processing unit is configured to perform cross-correlation calculations on the forward continuation and the residual backward continuation to obtain the initial multi-parameter imaging;

[0112] The fifth processing unit is configured to obtain OBS4C data from the initial multi-parameter imaging based on the Born forward modeling of the ACE equation, perform reverse data residual extrapolation on the simulated OBS4C data and the OBS4C data from the observation system, and perform cross-correlation calculation on the forward extrapolation and the reverse data residual extrapolation. After multiple iterations, the LSM imaging results are obtained, and then the seismic imaging results are obtained.

[0113] Example 3: This example provides an electronic device corresponding to the new elastic wave least squares reverse time migration method provided in Example 1. The electronic device can be an electronic device for the client, such as a mobile phone, laptop, tablet computer, desktop computer, etc., to execute the method of Example 1.

[0114] like Figure 7 As shown, the electronic device includes a processor, memory, a communication interface, and a bus. The processor, memory, and communication interface are connected via the bus to enable communication between them. The bus can be an Industry Standard Architecture (ISA) bus, a Peripheral Component Interconnect (PCI) bus, or an Extended Industry Standard Component (EISA) bus, etc. The memory stores a computer program that can run on the processor. When the processor runs the computer program, it executes the above-described method. The implementation principle and technical effects are similar to those in Embodiment 1, and will not be repeated here. Those skilled in the art will understand that... Figure 7 The structure shown is merely a block diagram of a portion of the structure related to the present application and does not constitute a limitation on the computing device on which the present application is applied. The specific computing device may include more or fewer components than shown in the figure, or combine certain components, or have different component arrangements.

[0115] In a preferred embodiment, the logical instructions in the aforementioned memory can be implemented as software functional units and sold or used as independent products, and can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, or a part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods described in the various embodiments of this application. The aforementioned storage medium includes various media capable of storing program code, such as USB flash drives, portable hard drives, read-only memory (ROM), random access memory (RAM), and optical discs.

[0116] In a preferred embodiment, the processor can be any type of general-purpose processor such as a central processing unit (CPU) or a digital signal processor (DSP), and is not limited thereto.

[0117] Example 4: This example provides a computer program product. The computer program product may include a computer program stored on a computer-readable storage medium. The computer program includes program instructions. When the program instructions are executed by the computer, the computer can execute the method provided in Example 1 above. Its implementation principle and technical effects are similar to those in Example 1, and will not be repeated here.

[0118] In a preferred embodiment, the computer-readable storage medium may be a tangible device for holding and storing instructions used by an instruction execution device, such as, but not limited to, an electrical storage device, a magnetic storage device, an optical storage device, an electromagnetic storage device, a semiconductor storage device, or any combination thereof. The computer-readable storage medium stores computer program instructions that cause a computer to perform the method provided in Embodiment 1 above.

[0119] The various embodiments in this specification are described in a progressive manner. Similar or identical parts between embodiments can be referred to interchangeably. Each embodiment focuses on describing the differences from other embodiments. In the description of this specification, the reference to terms such as "a preferred embodiment" indicates that a specific feature, structure, material, or characteristic described in connection with that embodiment or example is included in at least one embodiment or example of the embodiments in this specification. In this specification, the illustrative expressions of the above terms do not necessarily refer to the same embodiment or example. Furthermore, the specific features, structures, materials, or characteristics described can be combined in any suitable manner in one or more embodiments or examples. Moreover, without contradiction, those skilled in the art can combine and integrate the different embodiments or examples described in this specification, as well as the features of different embodiments or examples.

[0120] This application is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of this application. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart... Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0121] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0122] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0123] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A novel least-squares reverse-time migration method for elastic waves, characterized in that... include: Acquire P-wave velocity model, S-wave velocity model, and density model, as well as data from the OBS 4C observation system; Forward wavefield extension based on P-wave velocity model, S-wave velocity model and density model; Residual inverse extrapolation is performed based on the P-wave velocity model, S-wave velocity model, density model, and 4C data; The initial multi-parameter imaging is obtained by cross-correlation calculation of forward continuation and residual backward continuation; Initial multi-parameter imaging was performed using Born forward modeling based on the ACE equations to obtain OBS 4C data. The simulated OBS 4C data was then subjected to residual backward continuation with the observation system's OBS 4C data. Cross-correlation calculations were performed on the forward continuation and residual backward continuation, and after multiple iterations, LSM imaging results were obtained, thus yielding seismic imaging results. The formula used to obtain OBS 4C data from the initial multi-parameter imaging based on Born forward modeling based on the ACE equations is as follows: Where δ represents the disturbance, ( x , z (ξ, η) represents the Cartesian coordinates of the grid point before the coordinate transformation, and (ξ, η) represents the body-fitted coordinates of the grid point after the coordinate transformation. P It is acoustic stress. τ is density, t is time, and τ is time. xx and τ zz τ represents the normal stress component. xz Represents the shear stress component. v x0 , v z0 Representing the horizontal direction respectively x Direction, vertical direction z velocity components in the direction, and They are respectively and The disturbance λ , μ λ0 and μ0 are Lamé constants. v x , v z Representing the horizontal direction respectively x Direction, vertical direction z velocity components in the direction, , , and It is the directional partial derivative.

2. The novel elastic wave least-squares reverse-time migration method according to claim 1, characterized in that, The formula used for forward wavefield extension based on the P-wave velocity model, S-wave velocity model, and density model is as follows: in, f For the focal term, τ xx0 and τ zz0 τ represents the normal stress component. xz0 P represents the shear stress component, and P0 is the acoustic stress.

3. The novel elastic wave least-squares reverse-time migration method according to claim 1, characterized in that, The formula for calculating the residual inverse extrapolation based on the P-wave velocity mode, S-wave velocity model, density model, and 4C data is as follows: Wherein, the superscript R represents the wave field value at the receiver location, x r Indicates the location of the detector point. and These represent the composite data in the horizontal and vertical directions, respectively. With observation data The difference between them.

4. The elastic wave least-squares reverse-time migration method according to claim 1, characterized in that, The LSM imaging results are obtained by cross-correlation calculation of forward continuation and residual backward continuation multiple times, and then the seismic imaging results are output, including: Construct an objective function and determine whether the objective function meets the set requirements. If it does, output the final seismic imaging result. Otherwise, apply the diagonal Hessian operator to update the LSM imaging result. Calculate the step size using the conjugate gradient method and perform multiple iterations until the set requirements are met, then output the final seismic imaging result.

5. The elastic wave least-squares reverse-time migration method according to claim 4, characterized in that, The gradient formula used to update LSM imaging results is: Where g represents the gradient direction, λ and μ τ is Lamé's constant, where the superscript R indicates the wavefield value at the receiver location. xx and τ zz Let E represent the normal stress component and E represent the objective function.

6. A system for implementing the novel elastic wave least-squares reverse-time migration method according to any one of claims 1-5, characterized in that, The system includes: The first processing unit is configured to acquire P-wave velocity model, S-wave velocity model and density model, as well as OBS 4C observation system data; The second processing unit is configured to perform forward wavefield extension based on the longitudinal wave velocity model, the transverse wave velocity model, and the density model. The third processing unit is configured to perform residual inverse extrapolation based on the P-wave velocity model, S-wave velocity model, density model, and 4C data. The fourth processing unit is configured to perform cross-correlation calculations on the forward continuation and the residual backward continuation to obtain the initial multi-parameter imaging; The fifth processing unit is configured to obtain OBS 4C data from the initial multi-parameter imaging based on the Born forward modeling of the ACE equation, perform reverse data residual extrapolation on the simulated OBS 4C data and the OBS 4C data from the observation system, and perform cross-correlation calculation on the forward extrapolation and the reverse data residual extrapolation. After multiple iterations, the LSM imaging results are obtained, and then the seismic imaging results are obtained.

7. A computer-readable storage medium for storing one or more programs, characterized in that, The one or more programs include instructions that, when executed by a computing device, cause the computing device to perform any of the methods described in claims 1 to 5.

8. An electronic device, characterized in that, include: One or more processors, a memory, and one or more programs, wherein the one or more programs are stored in the memory and configured to be executed by the one or more processors, the one or more programs including instructions for performing any of the methods described in claims 1 to 5.

Images