Anisotropic medium pure wave least squares reverse time migration method
By employing the pure acoustic wave least squares reverse time migration method for anisotropic media, the problems of low imaging accuracy and high computational load in existing technologies have been solved, enabling efficient and accurate seismic exploration imaging under complex geological conditions.
Patent Information
- Application Number
- CN202511172914.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-08-21
- Publication Date
- 2025-11-21
- Estimated Expiration
- 2045-08-21
AI Technical Summary
Existing least-squares reverse time migration methods cannot effectively handle the seismic wave characteristics in anisotropic media, resulting in low imaging accuracy, high computational cost, and poor stability, making them unsuitable for oil and gas exploration under complex geological conditions.
The least squares reverse time migration method for pure acoustic waves in anisotropic media is adopted. By calculating the inverse migration operator of the pure acoustic wave equation in anisotropic media, the adjoint equation and parameter perturbation gradient are determined. The Hessian matrix and vector product are solved, and the velocity and anisotropic parameter perturbation model are iteratively updated. The Gaussian-Newton gradient preconditioning is used to handle the model parameter perturbation.
It improves the imaging accuracy and convergence speed of seismic data in anisotropic media, provides efficient imaging analysis tools, and enhances the seismic exploration results under complex geological conditions.
Smart Images

Figure CN120686337B_ABST
Abstract
Description
Technical Field
[0001] This application relates to the fields of data processing and oil and gas geophysical exploration technology, and in particular to a pure acoustic least squares reverse time migration method for anisotropic media. Background Technology
[0002] As exploration progresses, the targets of oil and gas exploration become increasingly complex, posing new challenges to existing seismic exploration methods. Migration imaging is a crucial technical aspect of seismic exploration observation system design and seismic data processing. Optimizing acquisition parameters through imaging analysis can improve the exploration results of target layers. Simultaneously, migration imaging is also an effective means of acquiring subsurface structures and geological formations. The reverse-time migration method based on the two-way wave equation has no dip angle limitations and can image any steeply dipping interface, achieving higher accuracy than Koschkhov integral migration and one-way wave migration. However, reverse-time migration uses the adjoint approximation of the inverse operator of the forward modeling operator, and when seismic data is incomplete (e.g., noisy, band-limited, missing, etc.), satisfactory imaging results cannot be obtained. Existing technologies provide least-squares reverse-time migration as an inversion method for linear waveforms. By iteratively updating the reflection coefficient / model parameter perturbation profile through inversion, the quality of seismic data migration imaging can be further improved.
[0003] Anisotropic media are media whose physical properties (such as absorbance, refractive index, conductivity, tensile strength, etc.) vary with direction, unlike isotropic media. They are commonly found in materials such as crystals and wood. In the field of exploration, anisotropy is prevalent in subsurface media. The directional arrangement of subsurface rocks, minerals, and fractures causes changes in the propagation speed of seismic waves along different directions (resulting in anisotropy). The isotropic medium assumption is increasingly unable to meet the needs of actual production in exploration in complex areas. Therefore, achieving high-precision seismic wave imaging in anisotropic media has become crucial for the current exploration of deep, deep-water, and unconventional oil and gas reservoirs.
[0004] Furthermore, most of the existing least-squares reverse time migration methods mentioned above have the following problems:
[0005] (1) Based on the isotropic assumption, the anisotropic characteristics of seismic waves cannot be considered;
[0006] (2) Based on the assumption of elastic anisotropy, it is difficult to handle the coupling problem of longitudinal waves and fast and slow transverse waves, and the computational load is large;
[0007] (3) Based on the pseudo-acoustic anisotropy assumption, transverse wave artifact interference cannot be completely eliminated, and the stability is poor. However, pure acoustic anisotropic least squares reverse time migration can better overcome the above problems. The computational load and speed of the method are moderate, and it is well applicable to the design of acquisition and observation systems, with broad application prospects. Compared with isotropic media, anisotropic media involve more model parameters (such as velocity and Thomsen anisotropy parameters). Crosstalk between perturbations of different parameters seriously affects the accuracy and convergence of reverse time migration. Summary of the Invention
[0008] This application discloses a least-squares reverse time migration method for pure acoustic waves in anisotropic media.
[0009] In a first aspect, this application discloses a least-squares reverse-time migration method for pure acoustic waves in anisotropic media, the method comprising:
[0010] Calculate the inverse migration operator for the pure acoustic wave equation in anisotropic media;
[0011] Determine the adjoint equation and parameter perturbation gradient formula for the pure acoustic wave equation in anisotropic media;
[0012] Solve the pure acoustic wave equations for anisotropic media using the inverse migration operator and the adjoint equations;
[0013] Determine the Hessian matrix and vector product formula for the pure acoustic wave equation in anisotropic media;
[0014] Calculate the perturbation gradient of the model parameters and the preconditions for the Gaussian-Newton gradient;
[0015] Iterative update speed and anisotropic parameter perturbation model.
[0016] Optionally, the inverse migration operator of the pure acoustic equation for anisotropic media is calculated in the following manner:
[0017] The equations for pure acoustic waves in VTI vertically and laterally isotropic and anisotropic media are as follows:
[0018] (1)
[0019] in, For wave field, For model parameters, For the focal term;
[0020] (2)
[0021] in, Let be the propagation velocity along the axis of symmetry. and For Thomsen anisotropy parameters, For time, and For spatial coordinates, Represents the transpose of a matrix or vector;
[0022] Based on the Born approximation, the seismic wavefield and model parameters are decomposed into background components and disturbance components, and formula (1) becomes:
[0023] (3)
[0024] in, For the background wave field, To disturb the wave field, For background parameters, These are the disturbance parameters;
[0025] Expanding formula (3) using Taylor series yields:
[0026] (4)
[0027] in,
[0028] (5)
[0029] (6)
[0030] (7)
[0031] The background wave field and background parameters also satisfy the wave equation shown in formula (1):
[0032] (8)
[0033] Subtracting formula (8) from formula (3) and ignoring higher-order terms of wave field and parameter perturbations, we get:
[0034] (9)
[0035] To eliminate the dimension of velocity, the image of velocity is represented by the relative change of velocity, as shown in expression (10):
[0036] (10)
[0037] Formula (10) can be rewritten as:
[0038] (11)
[0039] Among them, formulas (8) and (11) are the inverse migration operators of the pure acoustic wave equation for anisotropic media in VTI, that is, the perturbation wave field is calculated by perturbing the known model parameters. .
[0040] Optionally, the adjoint equation and parameter perturbation gradient formula of the pure acoustic wave equation for anisotropic media are determined in the following manner:
[0041] Establish the minimum objective function:
[0042] (12)
[0043] in, To simulate a disturbed wave field, To observe the perturbation wave field, For the calculation area, To the maximum recording time, For the detector point projection operator;
[0044] Solving this constrained optimization problem using the Lagrange multiplier method, the objective function becomes:
[0045] (13)
[0046] in, For the accompanying wave field;
[0047] Integrating equation (13) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain:
[0048] (14)
[0049] make The adjoint equation of the pure acoustic wave equation for VTI anisotropic media can be obtained:
[0050] (15)
[0051] Based on the chain rule Derive the gradient of the objective function with respect to the perturbation of the model parameters:
[0052] (16)
[0053] (17)
[0054] (18)
[0055] in, .
[0056] Optionally, the inverse migration operator and adjoint equation of the pure acoustic wave equation for anisotropic media are solved in the following manner:
[0057] The hyperbolic differential equations in formulas (8), (11) and (15) are solved quickly using high-order regular grid finite difference. The mixed partial derivatives in the x and z directions are approximated by finite difference along the two directions respectively. The constant coefficient Poisson equations in formulas (8), (11) and (15) are solved efficiently using a fast Poisson equation solver.
[0058] Optionally, the formula for the Hessian matrix and vector product of the pure acoustic wave equation for anisotropic media is determined in the following manner:
[0059] Model parameter perturbation in discrete cases The gradient of the objective function with respect to the perturbation of the model parameters and Hessian matrix They are represented as follows:
[0060] (19)
[0061] (20)
[0062] (twenty one)
[0063] (twenty two)
[0064] Where N is the dimension of the discrete grid;
[0065] gradient The first derivative of the objective function with respect to the perturbation of the model parameters and the Hessian matrix are represented. Let represent the second derivative of the objective function with respect to the perturbation of the model parameters. The relationship between the two is shown in expression (23):
[0066] (twenty three)
[0067] Construct a new objective function F as shown in expression (24):
[0068] (twenty four)
[0069] Where x is an arbitrary column vector of dimension 3N.
[0070] (25)
[0071] Based on formulas (23) and (24), we can obtain:
[0072] (26)
[0073] In the continuous case, the derivative formula (24) of the objective function with respect to the perturbation of the model parameters becomes:
[0074] (27)
[0075] in, ( () is any function related to spatial location;
[0076] Based on the Lagrange multiplier method, the derivative formula of the objective function F with respect to the perturbation of the model parameters is derived, and formula (27) becomes:
[0077] (28)
[0078] in, ( ), and These are Lagrange multiplier functions;
[0079] Formula (28) can be simplified to:
[0080] (29)
[0081] make , (30)
[0082] Formula (29) degenerates into:
[0083] (31)
[0084] Integrating equation (31) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain:
[0085] (32)
[0086] make and We can obtain:
[0087] (33)
[0088] (34)
[0089] Based on the chain rule Introduce the gradient of the new objective function with respect to the perturbation of the model parameters:
[0090] (35)
[0091] (36)
[0092] (37)
[0093] in, In the discrete case, the gradient of F with respect to the perturbation of the model parameters m is the product of the Hessian matrix H and any vector x.
[0094] Optionally, the perturbation gradients of the model parameters and the Gaussian-Newton gradient preconditions are calculated in the following manner, including:
[0095] The background wave field is obtained by numerically solving the inverse migration operator along the forward time direction according to formulas (8) and (11). and perturbation wave field The adjoint wave field is obtained by numerically solving the adjoint equation in reverse time according to formula (15). Based on formulas (16) to (18), calculate the gradient of the objective function with respect to the perturbation of the model parameters. ( );
[0096] Gradient of model parameter perturbation using the Hessian operator Perform preconditioning:
[0097] (38)
[0098] in, It is the inverse of the Hessian matrix. The gradient of the perturbation parameters of the model after preconditioning;
[0099] Directly calculating and storing the Hessian matrix or its inverse is quite difficult. Therefore, the inversion operation in formula (38) is transformed into solving the following system of linear equations:
[0100] (39)
[0101] The steps for solving the system of equations shown in formula (38) using the conjugate gradient method are as follows:
[0102] (1) Given initial values , , and ;
[0103] (2) Calculate the Hessian matrix using formulas (8), (30), (33), (34), (35)-(37). With vector product ;
[0104] (3) Calculation and : and ;
[0105] (4) Update and : , , , ;
[0106] (5) Repeat steps (2) to (4) until the set maximum number of iterations is reached, and output the final gradient of the model perturbation parameters. .
[0107] Optionally, the velocity and anisotropic parameter perturbation model is iteratively updated in the following manner:
[0108] gradient after preconditioning The velocity and anisotropic parameter perturbation model is updated, and the iteration steps are repeated until the convergence condition is met, resulting in the updated formula as follows:
[0109] (40)
[0110] in, For the number of iterations, The model parameters are perturbed for the current iteration and the next iteration. The model parameters are perturbed for the next iteration. This is the iteration step size.
[0111] In a second aspect, this application discloses an electronic device comprising: a processor, a memory, and a computer program stored in the memory and executable on the processor, wherein the computer program, when executed by the processor, implements the method as described in any of the preceding claims.
[0112] Thirdly, this application discloses a computer-readable storage medium on which a computer program is stored, which, when executed by a processor, implements the method as described in any of the preceding claims.
[0113] Fourthly, this application discloses a computer program product in which, when the instructions in the computer program product are executed by a processor of an electronic device, the electronic device implements the method described in any of the preceding claims.
[0114] The technical solution provided in this application may include the following beneficial effects:
[0115] By employing a novel truncated Gaussian-Newton Hessian operator gradient preconditioning to suppress crosstalk between perturbations of different model parameters, a reliable velocity and anisotropic parameter perturbation model is established. This provides an efficient and accurate imaging analysis tool for seismic exploration acquisition design in complex areas, improves the migration imaging effect of seismic data in anisotropic media, and enhances the imaging accuracy and convergence speed of pure acoustic least-squares reverse-time migration in anisotropic media. Furthermore, this application can obtain accurate velocity and anisotropic parameter perturbation models, which helps to improve the quality of seismic data migration imaging under complex geological conditions. Attached Figure Description
[0116] Figure 1 A flowchart of a method for the least-squares reverse-time migration of pure acoustic waves in anisotropic media provided in this application;
[0117] Figure 2 A schematic diagram of the Hess VTI model provided in this application under the velocity v parameter;
[0118] Figure 3 The Hess VTI model provided in this application has anisotropic parameters The diagram below;
[0119] Figure 4 The Hess VTI model provided in this application has anisotropic parameters The diagram below;
[0120] Figure 5 Velocity perturbation of the Hess VTI model using conventional least squares reverse time migration results A schematic diagram;
[0121] Figure 6 Anisotropic parameter perturbation of the Hess VTI model using conventional least-squares reverse time migration results A schematic diagram;
[0122] Figure 7 Anisotropic parameter perturbation of the Hess VTI model using conventional least-squares reverse time migration results A schematic diagram;
[0123] Figure 8 The velocity perturbation of the Hess VTI model using the least squares reverse time migration results provided in this application. A schematic diagram;
[0124] Figure 9 Anisotropic parameter perturbation of the Hess VTI model using the least squares reverse time migration results provided in this application. A schematic diagram;
[0125] Figure 10 Anisotropic parameter perturbation of the Hess VTI model using the least squares reverse time migration results provided in this application. A schematic diagram;
[0126] Figure 11 A block diagram of an electronic device provided in this application;
[0127] Figure 12 A block diagram of another electronic device provided in this application. Detailed Implementation
[0128] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, not all embodiments. Based on the embodiments of this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.
[0129] This application belongs to the fields of data processing and oil and gas geophysical exploration technology, and relates to a pure acoustic least squares reverse time migration method for anisotropic media. It can be applied to seismic acquisition design and high-precision imaging of seismic data under complex geological conditions, and can acquire velocity and anisotropic parameter perturbation profiles. The pure acoustic least squares reverse time migration method for anisotropic media provided in this application is described in detail below. Example 1
[0130] Reference Figure 1 This is a flowchart of a method for least-squares reverse-time migration of pure acoustic waves in anisotropic media, provided in this application. The method includes:
[0131] Step S101: Calculate the inverse migration operator for the pure acoustic wave equation of anisotropic media;
[0132] Step S102: Determine the adjoint equation and parameter perturbation gradient formula for the pure acoustic wave equation of the anisotropic medium;
[0133] Step S103: Solve the pure acoustic wave equation, inverse migration operator, and adjoint equation for the anisotropic medium;
[0134] Step S104: Determine the Hessian matrix and vector product formula for the pure acoustic wave equation of anisotropic media;
[0135] Step S105: Calculate the perturbation gradient of the model parameters and the Gaussian-Newton gradient preconditions;
[0136] Step S106: Iteratively update the perturbation model of the velocity and anisotropic parameters.
[0137] This application employs a novel truncated Gaussian-Newton Hessian operator gradient preconditioning to suppress crosstalk between perturbations of different model parameters, establishing a reliable velocity and anisotropic parameter perturbation model. This provides an efficient and accurate imaging analysis tool for seismic exploration acquisition design in complex areas, improves the migration imaging effect of seismic data in anisotropic media, and enhances the imaging accuracy and convergence speed of pure acoustic least-squares reverse-time migration in anisotropic media. Furthermore, this application can obtain accurate velocity and anisotropic parameter perturbation models, contributing to the improvement of the quality of seismic data migration imaging under complex geological conditions.
[0138] For clarity, the specific process of each step in Example 1 will be explained separately.
[0139] In step S101, the inverse migration operator of the pure acoustic wave equation for anisotropic media is calculated as follows:
[0140] The equations for pure acoustic waves in VTI vertically and laterally isotropic and anisotropic media are as follows:
[0141] (1)
[0142] in, For wave field, For model parameters, For the focal term;
[0143] (2)
[0144] Please see below. Figures 2 to 4 , Let be the propagation velocity along the axis of symmetry. and For Thomsen anisotropy parameters, For time, and For spatial coordinates, Represents the transpose of a matrix or vector;
[0145] Based on the Born approximation, the seismic wavefield and model parameters are decomposed into background components and disturbance components, and formula (1) becomes:
[0146] (3)
[0147] in, For the background wave field, To disturb the wave field, For background parameters, These are the disturbance parameters;
[0148] Expanding formula (3) using Taylor series yields:
[0149] (4)
[0150] in,
[0151] (5)
[0152] (6)
[0153] (7)
[0154] The background wave field and background parameters also satisfy the wave equation shown in formula (1):
[0155] (8)
[0156] Subtracting formula (8) from formula (3) and ignoring higher-order terms of wave field and parameter perturbations, we get:
[0157] (9)
[0158] To eliminate the dimension of velocity, the image of velocity is represented by the relative change of velocity, as shown in expression (10):
[0159] (10)
[0160] Formula (10) can be rewritten as:
[0161] (11)
[0162] Among them, formulas (8) and (11) are the inverse migration operators of the pure acoustic wave equation for anisotropic media in VTI, that is, the perturbation wave field is calculated by perturbing the known model parameters. It should be noted that, under the Born approximation, the background wave field... and They can be viewed as direct waves and primary reflected waves, respectively.
[0163] In step S102, the adjoint equation and parameter perturbation gradient formula of the pure acoustic wave equation for the anisotropic medium are determined in the following manner:
[0164] Establish the minimum objective function:
[0165] (12)
[0166] in, To simulate a disturbed wave field, To observe the perturbation wave field, For the calculation area, To the maximum recording time, For the detector point projection operator;
[0167] Solving this constrained optimization problem using the Lagrange multiplier method, the objective function becomes:
[0168] (13)
[0169] in, For the accompanying wave field;
[0170] Integrating equation (13) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain:
[0171] (14)
[0172] make The adjoint equation of the pure acoustic wave equation for VTI anisotropic media can be obtained:
[0173] (15)
[0174] Based on the chain rule Derive the gradient of the objective function with respect to the perturbation of the model parameters:
[0175] (16)
[0176] (17)
[0177] (18)
[0178] in, .
[0179] In step S103, the inverse migration operator and adjoint equation of the pure acoustic wave equation for the anisotropic medium are solved in the following manner:
[0180] The hyperbolic differential equations in formulas (8), (11) and (15) are solved quickly using high-order regular grid finite difference. The mixed partial derivatives in the x and z directions are approximated by finite difference along the two directions respectively. The constant coefficient Poisson equations in formulas (8), (11) and (15) are solved efficiently using a fast Poisson equation solver, which has higher computational efficiency than the pseudospectral method.
[0181] In step S104, the Hessian matrix and vector product formula for the pure acoustic wave equation of the anisotropic medium are determined in the following manner:
[0182] Model parameter perturbation in discrete cases The gradient of the objective function with respect to the perturbation of the model parameters and Hessian matrix They are represented as follows:
[0183] (19)
[0184] (20)
[0185] (twenty one)
[0186] (twenty two)
[0187] Where N is the dimension of the discrete grid;
[0188] gradient The first derivative of the objective function with respect to the perturbation of the model parameters and the Hessian matrix are represented. Let represent the second derivative of the objective function with respect to the perturbation of the model parameters. The relationship between the two is shown in expression (23):
[0189] (twenty three)
[0190] Construct a new objective function F as shown in expression (24):
[0191] (twenty four)
[0192] Where x is an arbitrary column vector of dimension 3N.
[0193] (25)
[0194] Based on formulas (23) and (24), we can obtain:
[0195] (26)
[0196] It should be noted that Equation (26) shows that the product of the Hessian matrix and any vector can be converted into the derivative of the objective function F with respect to the perturbation of the model parameters m.
[0197] In the continuous case, the derivative formula (24) of the objective function with respect to the perturbation of the model parameters becomes:
[0198] (27)
[0199] in, ( () is any function related to spatial location;
[0200] It should be noted that, due to the perturbation gradient of the model parameters... ( The simulated perturbation wave field is obtained through formulas (16)-(18). and accompanying wave field By solving the inverse offset operator shown in formula (11) and the adjoint equation shown in formula (5), formulas (11), (15) and (16) to (18) serve as constraints for the differentiation operation of F.
[0201] Based on the Lagrange multiplier method, the derivative formula of the objective function F with respect to the perturbation of the model parameters is derived, and formula (27) becomes:
[0202] (28)
[0203] in, ( ), and These are Lagrange multiplier functions;
[0204] Formula (28) can be simplified to:
[0205] (29)
[0206] make , (30)
[0207] Formula (29) degenerates into:
[0208] (31)
[0209] Integrating equation (31) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain:
[0210] (32)
[0211] make and We can obtain:
[0212] (33)
[0213] (34)
[0214] Based on the chain rule Introduce the gradient of the new objective function with respect to the perturbation of the model parameters:
[0215] (35)
[0216] (36)
[0217] (37)
[0218] in, In the discrete case, the gradient of F with respect to the perturbation of the model parameters m is the product of the Hessian matrix H and any vector x.
[0219] It should be noted that formulas (35) to (37) give the gradient formula of the new objective function F with respect to the perturbation of the model parameters in the continuous case. In the discrete case, it is only necessary to substitute the wave field and parameters of different grid points into the calculation. As can be seen from formula (26), the gradient of F with respect to the perturbation of the model parameters m in the discrete case is the product of the Hessian matrix H and any vector x.
[0220] In step S105, the perturbation gradient of the model parameters and the Gaussian-Newton gradient preconditions are calculated in the following manner, including:
[0221] The background wave field is obtained by numerically solving the inverse migration operator along the forward time direction according to formulas (8) and (11). and perturbation wave field The adjoint wave field is obtained by numerically solving the adjoint equation in reverse time according to formula (15). Based on formulas (16) to (18), calculate the gradient of the objective function with respect to the perturbation of the model parameters. ( );
[0222] Gradient of model parameter perturbation using the Hessian operator Perform preconditioning:
[0223] (38)
[0224] in, It is the inverse of the Hessian matrix. The gradient of the perturbation parameters of the model after preconditioning;
[0225] Directly calculating and storing the Hessian matrix or its inverse is quite difficult. Therefore, the inversion operation in formula (38) is transformed into solving the following system of linear equations:
[0226] (39)
[0227] The steps for solving the system of equations shown in formula (38) using the conjugate gradient method are as follows:
[0228] (1) Given initial values , , and ;
[0229] (2) Calculate the Hessian matrix using formulas (8), (30), (33), (34), (35)-(37). With vector product ;
[0230] (3) Calculation and : and ;
[0231] (4) Update and : , , , ;
[0232] (5) Repeat steps (2) to (4) until the set maximum number of iterations is reached, and output the final gradient of the model perturbation parameters. Preferably, the maximum number of iterations can be set to 5.
[0233] It should be noted that the above process does not require analytical calculation of the Hessian matrix or its inverse matrix. The Gaussian-Newton preconditioning of the perturbation parameter gradient is achieved by performing multiple multiplication operations (Hx) of the Hessian matrix and the vector.
[0234] In step S106, the velocity and anisotropic parameter perturbation model is iteratively updated in the following manner:
[0235] gradient after preconditioning The velocity and anisotropic parameter perturbation model is updated, and the iteration steps are repeated until the convergence condition is met (such as the value of the objective function being less than a certain threshold). The updated formula is as follows:
[0236] (40)
[0237] in, For the number of iterations, The model parameters are perturbed for the current iteration and the next iteration. The model parameters are perturbed for the next iteration. This is the iteration step size.
[0238] In one implementation, the iteration step size can be obtained through methods such as linear search or parabolic fitting.
[0239] This application proposes a novel least-squares reverse-time migration method for pure acoustic waves in anisotropic media. Based on the Born approximation, the inverse migration operator of the pure acoustic wave equation for anisotropic media is derived; the adjoint equation and parameter perturbation gradient formula of the pure acoustic wave equation for anisotropic media are derived using the Lagrange multiplier method; efficient solutions to the inverse migration operator and adjoint equation of the pure acoustic wave equation for anisotropic media are achieved using finite difference and fast Poisson solvers; the Hessian matrix and vector product formula of the pure acoustic wave equation for anisotropic media is derived based on the Lagrange multiplier method; the perturbation gradient of the model parameters is calculated and Gaussian-Newton gradient preprocessing is performed; and the velocity and anisotropic parameter perturbation model are iteratively updated.
[0240] This application employs a truncated Gauss-Newton inversion algorithm to mitigate the crosstalk effect caused by perturbations of different model parameters in anisotropic media. This helps to improve the convergence rate and imaging accuracy of pure acoustic least squares reverse time migration in anisotropic media, thereby improving the quality of oil and gas exploration under complex geological conditions.
[0241] Please see Figures 2 to 4 The Hess VTI model is used below to verify the effectiveness of the anisotropic medium pure acoustic wave least squares reverse time migration method proposed in this application.
[0242] Specifically, the model is divided into 630 × 218 grid points, with a grid size of 15.0m × 15.0m, a time step of 1.5ms, and a maximum recording time of 3.0s. The source uses a Ricker wavelet with a dominant frequency of 15.0 Hz. 63 shots and 630 receivers are evenly distributed on the surface, with a shot spacing of 150m and a trace spacing of 15m. Conventional anisotropic medium pure acoustic wave least squares reverse time migration method requires 8 forward models per iteration (2 forward extensions of the background wavefield and perturbation wavefield; 1 backward extension of the adjoint wavefield; 1 background wavefield reconstruction; and 4 iterations for calculating the iteration step size). The anisotropic medium pure acoustic wave least squares reverse time migration method of this application requires 28 forward models per iteration (8 for the conventional method; 4 for the gradient Hessian operator preconditioning, and 5 for the first iteration).
[0243] It should be noted that the conventional method requires the same amount of computation for 35 iterations as the method in this application requires 10 iterations (both are 280 forward modeling iterations). Figures 5 to 7 The results of 35 iterations of the conventional least-squares reverse time migration method are presented. As shown in the figure, the image with velocity perturbation is better than the image with anisotropic parameter perturbation.
[0244] Figures 8 to 10 The results of 10 iterations of the least-squares reverse-time migration method presented in this application are given. Low-frequency noise is effectively suppressed in the anisotropic parameter perturbation results, and the continuity of the in-phase axis is significantly improved. Compared with conventional methods, the velocity perturbation images obtained by the method in this application show better interface focusing and reduced ghosting.
[0245] The above results show that by adopting the new truncated Gauss-Newton inversion method, the pure acoustic least squares reverse time migration method for anisotropic media proposed in this application can provide a high-quality velocity and anisotropic parameter perturbation model, which significantly improves the accuracy and resolution of anisotropic media seismic data migration imaging.
[0246] It should be noted that, for the sake of simplicity, the method embodiments are all described as a series of actions. However, those skilled in the art should understand that this application is not limited to the described order of actions, because according to this application, some steps can be performed in other orders or simultaneously. Secondly, those skilled in the art should also understand that the embodiments described in the specification are all optional embodiments, and the actions involved are not necessarily required by this application.
[0247] Example 2
[0248] Optionally, this application also provides an electronic device, including: a processor, a memory, and a computer program stored in the memory and executable on the processor. When the computer program is executed by the processor, it implements the various processes of the above method embodiments and achieves the same technical effect. To avoid repetition, it will not be described again here.
[0249] This application also provides a computer-readable storage medium storing a computer program. When the computer program is executed by a processor, it implements the various processes of the above-described method embodiments and achieves the same technical effects. To avoid repetition, it will not be described again here. The computer-readable storage medium may be a read-only memory (ROM), a random access memory (RAM), a magnetic disk, or an optical disk, etc.
[0250] Figure 11 This application provides a block diagram of an electronic device 800. For example, the electronic device 800 may be a mobile phone, computer, digital broadcasting terminal, messaging device, game console, tablet device, medical device, fitness equipment, personal digital assistant, etc.
[0251] Reference Figure 11 The electronic device 800 may include one or more of the following components: a processing component 802, a memory 804, a power supply component 806, a multimedia component 808, an audio component 810, an input / output (I / O) interface 812, a sensor component 814, and a communication component 816.
[0252] Processing component 802 typically controls the overall operation of electronic device 800, such as operations associated with display, telephone calls, data communication, camera operation, and recording operations. Processing component 802 may include one or more processors 820 to execute instructions to complete all or part of the steps of the methods described above. Furthermore, processing component 802 may include one or more modules to facilitate interaction between processing component 802 and other components. For example, processing component 802 may include a multimedia module to facilitate interaction between multimedia component 808 and processing component 802.
[0253] Memory 804 is configured to store various types of data to support the operation of device 800. Examples of this data include instructions for any application or method operating on electronic device 800, contact data, phonebook data, messages, images, videos, etc. Memory 804 can be implemented by any type of volatile or non-volatile storage device or a combination thereof, such as static random access memory (SRAM), electrically erasable programmable read-only memory (EEPROM), erasable programmable read-only memory (EPROM), programmable read-only memory (PROM), read-only memory (ROM), magnetic storage, flash memory, magnetic disk, or optical disk.
[0254] Power supply component 806 provides power to various components of electronic device 800. Power supply component 806 may include a power management system, one or more power supplies, and other components associated with generating, managing, and distributing power to electronic device 800.
[0255] Multimedia component 808 includes a screen that provides an output interface between the electronic device 800 and the user. In some embodiments, the screen may include a liquid crystal display (LCD) and a touch panel (TP). If the screen includes a touch panel, the screen may be implemented as a touchscreen to receive input signals from the user. The touch panel includes one or more touch sensors to sense touches, swipes, and gestures on the touch panel. The touch sensors may sense not only the boundaries of the touch or swipe action but also the duration and pressure associated with the touch or swipe operation. In some embodiments, multimedia component 808 includes a front-facing camera and / or a rear-facing camera. When the device 800 is in an operating mode, such as a shooting mode or a video mode, the front-facing camera and / or the rear-facing camera may receive external multimedia data. Each front-facing camera and rear-facing camera may be a fixed optical lens system or have focal length and optical zoom capabilities.
[0256] Audio component 810 is configured to output and / or input audio signals. For example, audio component 810 includes a microphone (MIC) configured to receive external audio signals when electronic device 800 is in an operating mode, such as call mode, recording mode, and voice recognition mode. The received audio signals may be further stored in memory 804 or transmitted via communication component 816. In some embodiments, audio component 810 also includes a speaker for outputting audio signals.
[0257] I / O interface 812 provides an interface between processing component 802 and peripheral interface modules, such as keyboards, click wheels, buttons, etc. These buttons may include, but are not limited to, home buttons, volume buttons, power buttons, and lock buttons.
[0258] Sensor assembly 814 includes one or more sensors for providing state assessments of various aspects of electronic device 800. For example, sensor assembly 814 may detect the on / off state of device 800, the relative positioning of components such as the display and keypad of electronic device 800, changes in position of electronic device 800 or a component of electronic device 800, the presence or absence of user contact with electronic device 800, orientation or acceleration / deceleration of electronic device 800, and temperature changes of electronic device 800. Sensor assembly 814 may include a proximity sensor configured to detect the presence of nearby objects without any physical contact. Sensor assembly 814 may also include a light sensor, such as a CMOS or CCD image sensor, for use in imaging applications. In some embodiments, sensor assembly 814 may also include an accelerometer, gyroscope, magnetometer, pressure sensor, or temperature sensor.
[0259] Communication component 816 is configured to facilitate wired or wireless communication between electronic device 800 and other devices. Electronic device 800 can access wireless networks based on communication standards, such as WiFi, carrier networks (such as 2G, 3G, 4G, or 5G), or combinations thereof. In one exemplary embodiment, communication component 816 receives broadcast signals or broadcast operation information from an external broadcast management system via a broadcast channel. In one exemplary embodiment, communication component 816 also includes a near-field communication (NFC) module to facilitate short-range communication. For example, the NFC module may be implemented based on radio frequency identification (RFID) technology, Infrared Data Association (IrDA) technology, ultra-wideband (UWB) technology, Bluetooth (BT) technology, and other technologies.
[0260] In an exemplary embodiment, the electronic device 800 may be implemented by one or more application-specific integrated circuits (ASICs), digital signal processors (DSPs), digital signal processing devices (DSPDs), programmable logic devices (PLDs), field-programmable gate arrays (FPGAs), controllers, microcontrollers, microprocessors, or other electronic components to perform the methods described above.
[0261] In an exemplary embodiment, a non-transitory computer-readable storage medium including instructions is also provided, such as a memory 804 including instructions, which can be executed by a processor 820 of an electronic device 800 to perform the above-described method. For example, the non-transitory computer-readable storage medium may be a ROM, random access memory (RAM), CD-ROM, magnetic tape, floppy disk, and optical data storage device, etc.
[0262] Example 3
[0263] Figure 12 A block diagram of another electronic device 1900 provided for this application. For example, electronic device 1900 may be provided as a server.
[0264] Reference Figure 12 The electronic device 1900 includes a processing component 1922, which further includes one or more processors, and memory resources represented by memory 1932 for storing instructions, such as application programs, that can be executed by the processing component 1922. The application programs stored in memory 1932 may include one or more modules, each corresponding to a set of instructions. Furthermore, the processing component 1922 is configured to execute instructions to perform the methods described above.
[0265] Electronic device 1900 may also include a power supply component 1926 configured to perform power management of electronic device 1900, a wired or wireless network interface 1950 configured to connect electronic device 1900 to a network, and an input / output (I / O) interface 1958. Electronic device 1900 can operate on an operating system stored in memory 1932, such as Windows Server™, Mac OS X™, Unix™, Linux™, FreeBSD™, or similar.
[0266] Example 4
[0267] Fourthly, this application discloses a computer program product in which, when the instructions in the computer program product are executed by a processor of an electronic device, the electronic device is enabled to perform the method described in any of the preceding aspects.
[0268] It should be noted that, in this document, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Unless otherwise specified, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes that element.
[0269] Through the above description of the embodiments, those skilled in the art can clearly understand that the methods of the above embodiments can be implemented by means of software plus necessary general-purpose hardware platforms. Of course, they can also be implemented by hardware, but in many cases the former is a better implementation method. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product is stored in a storage medium (such as ROM / RAM, magnetic disk, optical disk) and includes several instructions to cause a terminal (which may be a mobile phone, computer, server, air conditioner, or network device, etc.) to execute the methods described in the various embodiments of this application.
[0270] The embodiments of this application have been described above with reference to the accompanying drawings. However, this application is not limited to the specific embodiments described above. The specific embodiments described above are merely illustrative and not restrictive. Those skilled in the art can make many other forms under the guidance of this application without departing from the spirit and scope of the claims, and all of these forms are within the protection scope of this application.
[0271] Those skilled in the art will recognize that the units and algorithm steps of the various examples described in conjunction with the embodiments disclosed in this application can be implemented in electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application.
[0272] Those skilled in the art will clearly understand that, for the sake of convenience and brevity, the specific working processes of the systems, devices, and units described above can be referred to the corresponding processes in the foregoing method embodiments, and will not be repeated here.
[0273] In the embodiments provided in this application, it should be understood that the disclosed apparatus and methods can be implemented in other ways. For example, the apparatus embodiments described above are merely illustrative. For instance, the division of units is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be through some interfaces; the indirect coupling or communication connection between apparatuses or units may be electrical, mechanical, or other forms.
[0274] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.
[0275] In addition, the functional units in the various embodiments of this application can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit.
[0276] If the aforementioned functions are implemented as software functional units and sold or used as independent products, they 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 portion 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, ROM, RAM, magnetic disks, or optical disks.
[0277] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in this application should be included within the scope of protection of this application. Therefore, the scope of protection of this application should be determined by the scope of the claims.
Claims
1. A least-squares reverse-time migration method for pure acoustic waves in anisotropic media, characterized in that, The method includes: Calculate the inverse migration operator for the pure acoustic wave equation in anisotropic media; Determine the adjoint equation and parameter perturbation gradient formula for the pure acoustic wave equation in anisotropic media; Solve the pure acoustic wave equations for anisotropic media using the inverse migration operator and the adjoint equations; Determine the Hessian matrix and vector product formula for the pure acoustic wave equation in anisotropic media; Calculate the perturbation gradient of the model parameters and the preconditions for the Gaussian-Newton gradient; Iterative update speed and anisotropic parameter perturbation model; The inverse migration operator for the pure acoustic wave equation of anisotropic media is calculated in the following manner: The equations for pure acoustic waves in VTI vertically and laterally isotropic and anisotropic media are as follows: (1) in, For wave field, For model parameters, For the focal term; (2) in, Let be the propagation speed along the axis of symmetry. and For Thomsen anisotropy parameters, For time, and For spatial coordinates, Represents the transpose of a matrix or vector; Based on the Born approximation, the seismic wavefield and model parameters are decomposed into background components and disturbance components, and formula (1) becomes: (3) in, For the background wave field, To disturb the wave field, For background parameters, These are the disturbance parameters; Expanding formula (3) using Taylor series yields: (4) in, (5) (6) (7) The background wave field and background parameters also satisfy the wave equation shown in formula (1): (8) Subtracting formula (8) from formula (3) and ignoring higher-order terms of wave field and parameter perturbations, we get: (9) To eliminate the dimension of velocity, the image of velocity is represented by the relative change of velocity, as shown in expression (10): (10) Formula (10) can be rewritten as: (11) Among them, formulas (8) and (11) are the inverse migration operators of the pure acoustic wave equation for anisotropic media in VTI, that is, the perturbation wave field is calculated by perturbing the known model parameters. .
2. The least-squares reverse-time migration method for pure acoustic waves in anisotropic media according to claim 1, characterized in that, The adjoint equation and parameter perturbation gradient formula for the pure acoustic wave equation of anisotropic media are determined as follows: Establish the minimum objective function: (12) in, To simulate a disturbed wave field, To observe the perturbation wave field, For the calculation area, To the maximum recording time, For the detector point projection operator; Solving this constrained optimization problem using the Lagrange multiplier method, the objective function becomes: (13) in, For the accompanying wave field; Integrating equation (13) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain: (14) make The adjoint equation of the pure acoustic wave equation for VTI anisotropic media can be obtained: (15) Based on the chain rule Derive the gradient of the objective function with respect to the perturbation of the model parameters: (16) (17) (18) in, .
3. The least-squares reverse-time migration method for pure acoustic waves in anisotropic media according to claim 2, characterized in that, The inverse migration operator and adjoint equations of the pure acoustic wave equations for anisotropic media are solved as follows: The hyperbolic differential equations in formulas (8), (11) and (15) are solved quickly using high-order regular grid finite difference. The mixed partial derivatives in the x and z directions are approximated by finite difference along the two directions respectively. The constant coefficient Poisson equations in formulas (8), (11) and (15) are solved efficiently using a fast Poisson equation solver.
4. The least-squares reverse-time migration method for pure acoustic waves in anisotropic media according to claim 3, characterized in that, The Hessian matrix and vector product formulas for the pure acoustic wave equations in anisotropic media are determined as follows: Model parameter perturbation in discrete cases The gradient of the objective function with respect to the perturbation of the model parameters and Hessian matrix They are represented as follows: (19) (20) (21) (22) Where N is the dimension of the discrete grid; gradient The first derivative of the objective function with respect to the perturbation of the model parameters and the Hessian matrix are represented. Let represent the second derivative of the objective function with respect to the perturbation of the model parameters. The relationship between the two is shown in expression (23): (23) Construct a new objective function F as shown in expression (24): (24) Where x is an arbitrary column vector of dimension 3N. (25) Based on formulas (23) and (24), we can obtain: (26) In the continuous case, the derivative formula (24) of the objective function with respect to the perturbation of the model parameters becomes: (27) in, ( () is any function related to spatial location; Based on the Lagrange multiplier method, the derivative formula of the objective function F with respect to the perturbation of the model parameters is derived, and formula (27) becomes: (28) in, ( ), and These are Lagrange multiplier functions; Formula (28) can be simplified to: (29) make , (30) Formula (29) degenerates into: (31) Integrating equation (31) by parts and assuming that the wave field at the initial, final, and boundary points is zero, we can obtain: (32) make and We can obtain: (33) (34) Based on the chain rule Introduce the gradient of the new objective function with respect to the perturbation of the model parameters: (35) (36) (37) in, In the discrete case, the gradient of F with respect to the perturbation of the model parameters m is the product of the Hessian matrix H and any vector x.
5. The least-squares reverse-time migration method for pure acoustic waves in anisotropic media according to claim 4, characterized in that, The perturbation gradients and Gaussian-Newton gradient preconditions for the model parameters are calculated as follows: The background wave field is obtained by numerically solving the inverse migration operator along the forward time direction according to formulas (8) and (11). and perturbation wave field The adjoint wave field is obtained by numerically solving the adjoint equation in reverse time according to formula (15). Based on formulas (16) to (18), calculate the gradient of the objective function with respect to the perturbation of the model parameters. ( ); Gradient of model parameter perturbation using the Hessian operator Perform preconditioning: (38) in, It is the inverse of the Hessian matrix. The gradient of the perturbation parameters of the model after preconditioning; Directly calculating and storing the Hessian matrix or its inverse is quite difficult. Therefore, the inversion operation in formula (38) is transformed into solving the following system of linear equations: (39) The steps for solving the system of equations shown in formula (38) using the conjugate gradient method are as follows: (1) Given initial values , , and ; (2) Calculate the Hessian matrix using formulas (8), (30), (33), (34), (35)-(37). With vector product ; (3) Calculation and : and ; (4) Update and : , , , ; (5) Repeat steps (2) to (4) until the set maximum number of iterations is reached, and output the final gradient of the model perturbation parameters. .
6. The method for least-squares reverse-time migration of pure acoustic waves in anisotropic media according to claim 5, characterized in that, The velocity and anisotropic parameter perturbation model is iteratively updated in the following manner, including: gradient after preconditioning The velocity and anisotropic parameter perturbation model is updated, and the iteration steps are repeated until the convergence condition is met, resulting in the updated formula as follows: (40) in, For the number of iterations, The model parameters are perturbed for the current iteration and the next iteration. The model parameters are perturbed for the next iteration. This is the iteration step size.
7. An electronic device, characterized in that, include: A processor, a memory, and a computer program stored in the memory and executable on the processor, wherein the computer program, when executed by the processor, implements the method as described in any one of claims 1 to 6.
8. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores a computer program that, when executed by a processor, implements the method as described in any one of claims 1 to 6.
9. A computer program product, characterized in that, When the instructions in the computer program product are executed by the processor of the electronic device, the electronic device implements the method as described in any one of claims 1 to 6.