Full waveform inversion method and device based on self-accompanying VTI equation, medium and equipment

Abstract

The invention relates to a full waveform inversion method and device based on a self-accompanying VTI equation, a medium and equipment, and the method comprises the steps: splitting a VTI acoustic wave equation and a VTI elastic wave equation, and obtaining the split VTI acoustic wave equation and VTI elastic wave equation; constructing auxiliary variables p and q, and substituting the auxiliary variables p and q into the split VTI acoustic wave equation and VTI elastic wave equation to obtain a self-accompanying VTI acoustic wave equation and a self-accompanying VTI elastic wave equation; solving the self-accompanying VTI acoustic wave equation and the self-accompanying VTI elastic wave equation, and calculating to obtain a forward transmission wave field u and prediction data dpre; using the objective function and the prediction data dpre to calculate a residual error dres; a self-accompanying VTI acoustic wave equation and a self-accompanying VTI elastic wave equation are solved, and an accompanying wave field u * is calculated through back transmission of a residual error dres; and on the basis of a gradient formula of the target function to model parameters, calculating a gradient # imgabs0 # by using the forward propagation wave field u and the accompanying wave field u *, calculating a step length alpha, updating the model m through the gradient # imgabs1 # and an iterative updating formula, and judging whether the target function is converged or not.

Description

Technical Field

[0001] The invention relates to a full waveform inversion method, device, medium and equipment based on a self-adjoint VTI equation, belonging to the technical field of geological exploration. Background Art

[0002] Full waveform inversion of vertically symmetric transversely isotropic (VTI) media can effectively model the underground media in the VTI anisotropic region. Full waveform inversion is usually performed on a good initial model through iterative updates of local optimization algorithms to achieve high-precision inversion of underground media. The iterative update of the local optimization algorithm requires the gradient of the objective function as the update direction, and the calculation of the gradient requires accurate numerical simulation of the forward wave field and the adjoint wave field.

[0003] At present, the VTI elastic wave equation can be used to simulate the forward wave field of the VTI elastic wave field. In order to simulate the VTI acoustic wave field, the shear wave velocity along the symmetry axis in the elastic wave equation is usually set to 0 (Duveneck et al., 2008) to obtain the VTI acoustic anisotropic wave equation, and accurately realize the forward wave field simulation of the VTI acoustic wave field ( Figure 1 (a)). In the process of realizing full waveform inversion of VTI anisotropic media, the adjoint state method is usually used to obtain the gradient of the objective function with respect to the model parameters. The adjoint state method requires the derivation of the adjoint wave equation and the adjoint wave field. However, the VTI elastic wave equation and the VTI acoustic wave equation that are currently commonly used are not self-symmetric. Therefore, there are obvious problems with the adjoint wave equation obtained directly by transposing the wave equation operator. First, although the elastic wave adjoint wave equation of VTI elastic media can accurately realize the numerical simulation of the adjoint wave field, the free boundary conditions of the adjoint wave equation of the direct transposed wave equation operator are difficult to implement effectively. The main reason is that the physical meaning of the adjoint variables of the adjoint wave equation after transposition is unclear. Therefore, in the current references, most scholars only consider the absorbing boundary in the wave field simulation and do not implement the free boundary (Ren and Liu, 2016). Secondly, there are two problems with the VTI acoustic wave adjoint state equation derived by direct transposition: 1. When the adjoint equation does not consider the influence of the model parameters changing with space, the wave field simulated by the adjoint wave equation will have simulation instability problems at the mutation position of the anisotropic parameters (such as Figure 1 (b)); 2. When the adjoint equation considers the influence of model parameters changing with space, the adjoint wave field will have the problem of weak amplitude value in deep parts ( Figure 1 (d)). Summary of the invention

[0004] In view of the above technical problems, the present invention provides a full waveform inversion method, device, medium and equipment based on the self-adjoint VTI equation. The adjoint equation derived by the method will not have the problems of unstable accompanying wave field simulation and weak deep amplitude value; secondly, the accompanying wave field of the method can accurately realize the surface free boundary condition; finally, the method can use the same set of codes to realize the numerical simulation of the forward wave field and the accompanying wave field, thereby reducing the complexity of programming.

[0005] To achieve the above object, the present invention adopts the following technical solutions:

[0006] A full waveform inversion method based on self-adjoint VTI equation, comprising:

[0007] The normal stress τ containing the anisotropy parameter ε in the VTI acoustic wave equation and the VTI elastic wave equation xx Perform splitting to obtain split VTI acoustic wave equation and VTI elastic wave equation;

[0008] Construct auxiliary variables p and q and bring them into the split VTI acoustic wave equation and VTI elastic wave equation to obtain the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation;

[0009] Solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation to calculate the forward wave field u and the predicted data d pre ;

[0010] Using the objective function and prediction data d pre Calculate the residual d res ;

[0011] Solve the self-adjoint VTI acoustic wave equation and the self-adjoint VTI elastic wave equation, and back-propagate the residual d res Calculate the accompanying wave field u * ;

[0012] Based on the gradient formula of the objective function to the model parameters, the forward wave field u and the accompanying wave field u are used. * Computing Gradients

[0013] Calculate the step size α and pass the gradient The model m is updated with the iterative update formula to determine whether the objective function converges. If the objective function converges, the result is output. If the objective function does not converge, the next iteration is performed until the objective function converges.

[0014] In the full waveform inversion method based on the self-adjoint VTI equation, preferably, the specific formula of the split VTI acoustic wave equation is as follows:

[0015]

[0016] Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represents the normal stress in the horizontal and vertical directions respectively; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent the density, P-wave velocity and anisotropy parameter, respectively.

[0017] In the full waveform inversion method based on the self-adjoint VTI equation, preferably, the specific formula of the split VTI elastic wave equation is as follows:

[0018]

[0019] Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represent the normal stress in the horizontal and vertical directions respectively; τ xz represents the shear stress in the shear direction; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter, respectively.

[0020] The full waveform inversion method based on the self-accompanied VTI equation is preferably a self-accompanied VTI acoustic wave equation having the following specific formula:

[0021]

[0022] Among them, p and q represent auxiliary variables; τ xx2 represents the partial normal stress in the horizontal direction after splitting; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent density, P-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

[0023] The full waveform inversion method based on the self-accompanied VTI equation, preferably, the specific formula of the self-accompanied VTI elastic wave equation is as follows:

[0024]

[0025] Among them, p and q represent auxiliary variables; τ xx2 represents the partial normal stress in the horizontal direction after splitting; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

[0026] In the full waveform inversion method based on the self-adjoint VTI equation, preferably, the gradient formula of the objective function to the model parameters is as follows:

[0027]

[0028] Where, “T” represents the transpose of the matrix; m represents the model parameter; u represents the forward wave field, A represents the forward operator; u * represents the accompanying wave field; x represents the spatial coordinate of a point underground.

[0029] In the full waveform inversion method based on the self-adjoint VTI equation, preferably, the formula of the objective function is as follows:

[0030]

[0031] Among them, m represents the model parameter, d pre represents the predicted data generated by model m, d obs Represents observation data collected in the field.

[0032] A second aspect of the present invention provides a full waveform inversion device based on a self-adjoint VTI equation, comprising:

[0033] The first processing unit is used to process the normal stress τ containing the anisotropy parameter ε in the VTI acoustic wave equation and the VTI elastic wave equation. xx Perform splitting to obtain split VTI acoustic wave equation and VTI elastic wave equation;

[0034] The second processing unit is used to construct auxiliary variables p and q and bring them into the split VTI acoustic wave equation and VTI elastic wave equation to obtain a self-accompanied VTI acoustic wave equation and a self-accompanied VTI elastic wave equation;

[0035] The third processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and calculate the forward wave field u and the predicted data d pre ;

[0036] The fourth processing unit is used to use the objective function and the prediction data d pre Calculate the residual d res ;

[0037] The fifth processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and back propagate the residual d res Calculate the accompanying wave field u * ;

[0038] The sixth processing unit is used to use the forward wave field u and the accompanying wave field u based on the gradient formula of the objective function to the model parameters * Computing Gradients

[0039] The seventh processing unit is used to calculate the step size α and pass the gradient The model m is updated with the iterative update formula to determine whether the objective function converges. If the objective function converges, the result is output. If the objective function does not converge, the next iteration is performed until the objective function converges.

[0040] A third aspect of the present invention provides a computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of any one of the above-mentioned full waveform inversion methods based on the self-adjoint VTI equations.

[0041] A fourth aspect of the present invention provides a computer device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the processor executes the computer program, the steps of any one of the above-mentioned full waveform inversion methods based on the self-adjoint VTI equation are implemented.

[0042] The present invention adopts the above technical solution, which has the following advantages:

[0043] 1. In the full waveform inversion based on the self-adjoint VTI equation, the present invention firstly adjusts the normal stress τ xx The system is split and then auxiliary variables p and q are constructed to realize the self-adjoint of the VTI equation. This technology can reduce the programming complexity of full waveform inversion of VTI media and solve the problems of unstable accompanying wavefield simulation and weak deep amplitude value in the prior art. In addition, the forward wavefield and the accompanying wavefield can adopt the same free boundary implementation scheme.

[0044] 2. The present invention uses the normal stress τ of the VTI acoustic wave equation xxThe self-adjointness of the VTI acoustic wave equation is realized by splitting and constructing auxiliary variables p and q. The full waveform inversion based on the self-adjoint VTI acoustic wave equation can effectively avoid the simulation instability of the accompanying wave field and the weak amplitude value in the deep.

[0045] 3. The present invention uses the normal stress τ of the VTI elastic wave equation xx The self-adjointness of the VTI elastic wave equation is realized by splitting and constructing auxiliary variables p and q. The full waveform inversion based on the self-adjoint VTI elastic wave equation can effectively avoid the problem that the adjoint wave field cannot effectively implement the free boundary. BRIEF DESCRIPTION OF THE DRAWINGS

[0046] Figure 1 The forward and companion wave fields corresponding to the Marmousi2 model at 1700 milliseconds provided in one embodiment of the present invention; wherein, Figures (a) and (c) are the forward wave fields calculated by the VTI acoustic wave equation (Formula (2)); Figure (b) is the companion wave field calculated by the companion equation without considering the spatial variation of model parameters (Formula (8)); Figure (d) is the companion wave field calculated by the companion equation considering the spatial variation of model parameters (Formula (9)); Figure (e) is the forward wave field calculated by the self-accompanying VTI acoustic wave equation (Formula (15)) proposed by the present invention, and Figure (f) is the companion wave field calculated by the self-accompanying VTI acoustic wave equation (Formula (15)) proposed by the present invention. The source location is 6.5 km, the depth is 10 m, and the surface is a free boundary;

[0047] Figure 2 are the real model and initial model of Marmousi2; among them, Figures (a), (c), (e), (g) and (i) are the real P-wave velocity model, S-wave velocity model, density model and anisotropic parameter model respectively, and Figures (b), (d), (f), (h) and (j) are the corresponding initial models respectively;

[0048] Figure 3 The VTI acoustic wave full waveform inversion results of the Marmousi2 model, wherein Figure (a) is the P-wave velocity inverted by the prior art; Figure (b) is the P-wave velocity inverted by the present invention;

[0049] Figure 4 The VTI elastic wave full waveform inversion results of the Marmousi2 model, wherein (a) and (c) are the P-wave velocity and S-wave velocity inverted by the prior art, respectively; (b) and (d) are the P-wave velocity and S-wave velocity inverted by the present invention, respectively;

[0050] Figure 5 It is a flow chart of the full waveform inversion equation based on the self-adjoint VTI equation. DETAILED DESCRIPTION

[0051] In order to make the purpose, technical solution and advantages of the present invention clearer, the technical solution of the present invention is described clearly and completely below. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary persons in the art without creative work are within the scope of protection of the present invention.

[0052] Unless otherwise defined, the technical or scientific terms used in the present invention shall have the usual meanings understood by persons with ordinary skills in the field to which the present invention belongs. The words "first", "second", "third", "fourth" and similar terms used in the present invention do not indicate any order, quantity or importance, but are only used to distinguish different components. Words such as "include" or "comprise" and the like mean that the elements or objects appearing before the word include the elements or objects listed after the word and their equivalents, without excluding other elements or objects. Words such as "connect" or "connected" and the like are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect.

[0053] For ease of description, spatially relative terms may be used herein to describe the relationship of one element or feature relative to another element or feature as shown in the figures, such as "inside", "outside", "inner side", "outer side", "below", "above", etc. Such spatially relative terms are intended to include different orientations of the device in use or operation in addition to the orientation depicted in the figures.

[0054] Existing full waveform inversion usually uses the least squares objective function to fit the surface observation records and simulated prediction records to invert the underground medium parameters. Its objective function can be expressed as:

[0055]

[0056] Among them, m represents the model parameter, d pre represents the predicted data generated by model m, d obs Represents observation data collected in the field.

[0057] For full waveform inversion of acoustic waves in VTI media, the objective function (1) is constrained by the VTI acoustic wave equation:

[0058]

[0059] Among them, τ xx and τ zz Represents the normal stress in the horizontal and vertical directions respectively; v x and v zRespectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent density, P-wave velocity and anisotropy parameters respectively; s represents the earthquake source, and is loaded in the form of a pure P-wave source. The above formula (2) gives the two-dimensional VTI acoustic wave equation, and there is also a three-dimensional VTI acoustic wave equation. The two are similar in form, and there is a variable τ along the y-axis yy and v y The specific form will not be described here.

[0060] For full waveform inversion of elastic waves in VTI media, the objective function (1) is constrained by the VTI elastic wave equation:

[0061]

[0062] Among them, ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter, respectively; f is a newly defined variable, which can be accessed through v p and v s Calculated; τ xx and τ zz Represent the normal stress in the horizontal and vertical directions respectively; τ xz represents the shear stress in the shear direction; v x and v z Represent the particle vibration velocity in the horizontal and vertical directions respectively; s represents the earthquake source, and is loaded in the form of a pure longitudinal wave source. The same as the acoustic wave equation, the above formula (3) gives the two-dimensional VTI elastic wave equation, and there is also a three-dimensional VTI elastic wave equation. The two are similar in form, and there is a variable τ related to the y-axis yy , τ xy , τ yz and v y The specific form will not be described here.

[0063] Formulas (2) and (3) can both be expressed in matrix form as:

[0064] Au(x,t)=s(x s ,t) (4)

[0065] Where u represents the forward wave field, which includes multiple components, which are specifically defined in formulas (2) and (3); x represents the spatial coordinates of a point underground; s represents the earthquake source, and x s represents the spatial coordinates of the earthquake source, and t represents the time of seismic wave propagation. A represents the forward operator. For the VTI acoustic wave equation shown in formula (2), its forward operator A can be specifically expressed as:

[0066]

[0067] For the VTI elastic wave equation shown in formula (3), its forward operator A can be specifically expressed as:

[0068]

[0069] The forward wave field in the full waveform inversion of VTI media can be numerically simulated by formula (4). The companion wave field requires the corresponding companion equation for numerical simulation. By deducing through the companion state method, the companion equation corresponding to formula (4) can be expressed in matrix form as follows:

[0070] A T u * (x,t)=d res (x r ,t) (7)

[0071] Among them, u * represents the accompanying wave field, including and Four components, A T represents the adjoint operator of the adjoint wavefield, which is the transpose of equations (5) and (6); d res represents the data residual as the accompanying source, x r Represents the spatial coordinates of the detector.

[0072] Formula (7) shows that without considering the spatial variation of model parameters, the adjoint equation of formula (2) can be expressed as:

[0073]

[0074] in, and represents the four components of the accompanying wave field; δp represents the residual data of the pressure as the accompanying source. Similarly, considering the model parameters varying with space, the accompanying equation of formula (2) can be expressed as:

[0075]

[0076] The main difference between formula (8) and (9) lies in whether the corresponding spatial derivatives of the model parameters are obtained during the solution of the adjoint equation.

[0077] Through formulas (1), (4) and (7), the gradient of the objective function with respect to the model parameters can be expressed as:

[0078]

[0079] Where, “T” represents the transpose of the matrix; m represents the model parameter; u represents the forward wave field, A represents the forward operator; u *represents the accompanying wave field; x represents the spatial coordinate of a point underground.

[0080] The gradient and local optimization algorithm shown in formula (10) can be used to iteratively update the initial model parameters, thereby realizing the inversion of the underground model parameters, that is, to obtain the underground model parameters, including the specific values ​​of the P-wave and S-wave velocities and anisotropy parameters. The specific iterative update formula is as follows:

[0081]

[0082] In the formula, m i represents the model parameters of the i-th iteration, m i+1 represents the updated gradient.

[0083] From the forward operators shown in formulas (5) and (6), it can be seen that operator A is not self-symmetric. The forward operator A and its adjoint operator A T are not equal, so different codes are needed to realize the numerical simulation of the forward wave field and the adjoint wave field respectively during the implementation of full waveform inversion. At the same time, the adjoint equation derived from the full waveform inversion based on the VTI acoustic wave equation (Formula 2) has two problems: 1. When A does not consider the spatial variation of model parameters, the adjoint equation will simulate instability at the location where the anisotropic parameters have a sudden change, such as Figure 1 (b) as shown; 2. When A T When considering the model parameters changing with space, the accompanying wave field will have the problem of weak amplitude values ​​in the deep part, such as Figure 1 (d) In addition, the adjoint equation derived from the full waveform inversion of the VTI elastic wave equation (Equation 3) is difficult to implement accurately for free boundaries.

[0084] In view of the above shortcomings of the existing full waveform inversion of VTI media, the present invention proposes a full waveform inversion method based on the self-adjoint VTI equation, which can accurately realize the stable and correct numerical simulation of the forward wave field and the adjoint wave field in the full waveform inversion of VTI media under the condition of free surface ( Figure 1 (e) and 1(f)). At the same time, due to the self-adjoint nature of the forward operator, i.e., self-symmetry, the forward operator is the same as the adjoint operator. The numerical simulation of the forward wave field and the adjoint wave field in programming only requires one set of code to execute, which can reduce the complexity of the code implementation.

[0085] like Figure 1 As shown, the full waveform inversion method based on the self-adjoint VTI equation provided by the present invention includes the following specific steps:

[0086] In the present invention, in order to realize the self-adjoint of the forward operator A, the normal stress τ containing the anisotropy parameter ε in equations (2) and (3) is firstly xxThe VTI acoustic wave equation shown in formula (2) can be expressed as follows after splitting:

[0087]

[0088] Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represents the normal stress in the horizontal and vertical directions respectively; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent the density, P-wave velocity and anisotropy parameter, respectively.

[0089] The VTI elastic wave equation shown in formula (3) can be expressed as follows after decomposition:

[0090]

[0091] Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represent the normal stress in the horizontal and vertical directions respectively; τ xz represents the shear stress in the shear direction; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter, respectively.

[0092] Then the forward operator A is made self-adjoint by constructing auxiliary variables p and q:

[0093]

[0094] and

[0095]

[0096] Substituting formula (14) and (15) into formula (12) yields the self-adjoint VTI wave equation:

[0097]

[0098] Among them, p and q represent auxiliary variables; τ xx2 represents the partial normal stress in the horizontal direction after splitting; v x and v zRespectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent density, P-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

[0099] Substituting formula (14) and (15) into formula (13) yields the self-adjoint VTI elastic wave equation:

[0100]

[0101] Among them, p and q represent auxiliary variables; τ xx2 represents the partial normal stress in the horizontal direction after splitting; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

[0102] The present invention also comprises the steps of:

[0103] The forward wave field u and the predicted data d are calculated by solving the self-accompanying VTI acoustic wave (Formula (16)) or elastic wave (Formula (17)) pre ;

[0104] Through the objective function (formula (1)) and the predicted data d pre Calculate the residual d res ;

[0105] The residual d is back-propagated by solving the self-accompanying VTI acoustic wave (Equation (16)) and elastic wave (Equation (17)) res Calculate the accompanying wave field u *

[0106] By using formula (10) to use the forward wave field u and the accompanying wave field u * Computing Gradients

[0107] Calculate the step size α and pass the gradient Update the model m with formula (11) and determine whether the objective function converges. If the objective function converges, output the result. If the objective function does not converge, perform the next iteration until the objective function converges.

[0108] The self-accompanying VTI acoustic and elastic wave equations shown in equations (16) and (17) have forward operators A equal to their transpose A T Therefore, the full waveform inversion based on formula (16) and formula (17) has the same forward operator A and adjoint operator A T In the implementation of full waveform inversion, the same set of codes can be used to realize the numerical simulation of the forward wavefield and the companion wavefield, and the companion wavefield can use the same free boundary implementation scheme as the forward wavefield.

[0109] The following describes in detail through specific examples how the present invention can effectively realize full waveform inversion of acoustic waves and elastic waves in VTI media.

[0110] The example used in the present invention is the Marmousi2 model. The real longitudinal wave velocity model is as follows: Figure 2 As shown in (a), the real shear wave velocity model is as follows Figure 2 As shown in (c), the actual density model is Figure 2 As shown in (e), the true anisotropic parameter ε model is Figure 2 As shown in (g), the true anisotropy parameter δ is Figure 2 (i) The initial velocity models required for full waveform inversion are as follows: Figure 1 (b), 1(d) and 1(f).

[0111] The present invention only shows the implementation process of full waveform inversion of acoustic waves in VTI media. Figure 2 The initial P-wave velocity model shown in (b) is updated, the density is generated using the Gardner formula, and the anisotropy parameters are not updated during the inversion process using the initial model. The P-wave velocity inversion result obtained by the prior art is as follows: Figure 3 (a) shows that the longitudinal wave velocity inversion result of the present invention is as follows: Figure 3 By comparison, it can be seen that the prior art has the problem of weak deep amplitude value of the accompanying wave field, and the deep updating effect of the final inversion result is not as good as that of the present invention, while the inversion result of the present invention is also closer to the real P-wave velocity model.

[0112] The present invention is used in the process of implementing the elastic wave full waveform inversion of VTI medium. Figure 2 (b) and Figure 2 The initial P-wave velocity and S-wave velocity models shown in (d) are updated, the density is generated using the Gardner formula, and the anisotropy parameters are not updated during the inversion process using the initial model. The P-wave velocity and S-wave velocity inversion results obtained by the prior art are shown in Figures 1 and 2. Figure 4 (a) and 4 (c), the inversion results of the longitudinal wave velocity and the shear wave velocity obtained by the present invention are as follows Figure 4(b) and 4(d). Figure 4 It can be seen from the comparison between (a) and 4(c) that the adjoint wave field of the prior art cannot effectively implement the free boundary, and the final inversion result contains a lot of noise. Figure 4 (b) and 4(d) can accurately implement the free boundary because of their accompanying wave fields, so the final inversion results are closer to the true velocity model than the existing technology.

[0113] In the full waveform inversion based on the self-adjoint VTI equation, the present invention firstly adjusts the normal stress τ xx The system is split and then auxiliary variables p and q are constructed to realize the self-adjoint of the VTI equation. This technology can reduce the programming complexity of full waveform inversion of VTI media and solve the problems of unstable accompanying wavefield simulation and weak deep amplitude value in the prior art. In addition, the forward wavefield and the accompanying wavefield can adopt the same free boundary implementation scheme.

[0114] The method of the present invention has the following beneficial effects: firstly, the adjoint equation derived by the method will not have the problems of unstable adjoint wavefield simulation and weak deep amplitude value; secondly, the adjoint wavefield of the method can accurately realize the free boundary; finally, the method can use the same set of codes to realize the numerical simulation of the forward wavefield and the adjoint wavefield, thereby reducing the complexity of programming.

[0115] In the technical solution of the present invention, the objective function of the full waveform inversion is an objective function based on least squares fitting. The objective function may also be other publicly published objective functions, such as the optimal transmission objective function, the normalized objective function, etc. However, these other objective functions still need to be combined with the self-adjoint equations in the present invention to achieve stable and efficient inversion. In addition, in order to simplify the derivation, the present invention only writes the derivation of the self-adjoint VTI acoustic wave equation and elastic wave equation in the two-dimensional case, and the self-adjoint can also be achieved through the above process for the three-dimensional case.

[0116] A second aspect of the present invention provides a full waveform inversion device based on a self-adjoint VTI equation, comprising:

[0117] The first processing unit is used to process the normal stress τ containing the anisotropy parameter ε in the VTI acoustic wave equation and the VTI elastic wave equation. xx Perform splitting to obtain split VTI acoustic wave equation and VTI elastic wave equation;

[0118] The second processing unit is used to construct auxiliary variables p and q and bring them into the split VTI acoustic wave equation and VTI elastic wave equation to obtain a self-accompanied VTI acoustic wave equation and a self-accompanied VTI elastic wave equation;

[0119] The third processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and calculate the forward wave field u and the predicted data d pre ;

[0120] The fourth processing unit is used to use the objective function and the prediction data d pre Calculate the residual d res ;

[0121] The fifth processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and back propagate the residual d res Calculate the accompanying wave field u * ;

[0122] The sixth processing unit is used to use the forward wave field u and the accompanying wave field u based on the gradient formula of the objective function to the model parameters * Computing Gradients

[0123] The seventh processing unit is used to calculate the step size α and pass the gradient The model m is updated with the iterative update formula to determine whether the objective function converges. If the objective function converges, the result is output. If the objective function does not converge, the next iteration is performed until the objective function converges.

[0124] A third aspect of the present invention provides a computer-readable storage medium having a computer program stored thereon, wherein the computer program, when executed by a processor, implements the steps of any one of the above-mentioned full waveform inversion methods based on the self-adjoint VTI equations.

[0125] A fourth aspect of the present invention provides a computer device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein when the processor executes the computer program, the steps of any one of the above-mentioned full waveform inversion methods based on the self-adjoint VTI equation are implemented.

[0126] The present invention is described in terms of flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to specific embodiments. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0127] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0128] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process in the computer or other programmable device. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0129] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some of the technical features therein. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A full waveform inversion method based on the self-adjoint VTI equation, characterized in that: include: The normal stress τ in the VTI acoustic wave equation and the VTI elastic wave equation xx Split into τ xx1 and τ xx2 , and obtain the split VTI acoustic wave equation and VTI elastic wave equation; Construct auxiliary variables p and q and bring them into the split VTI acoustic wave equation and VTI elastic wave equation to obtain the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation; Solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation to calculate the forward wave field u and the predicted data d pre ; Using the objective function and prediction data d pre Calculate the residual d res ; Solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and back propagate the residual d res Calculate the accompanying wave field u * ; Based on the gradient formula of the objective function to the model parameters, the forward wave field u and the accompanying wave field u are used. * Computing Gradients Calculate the step size α and pass the gradient The model m is updated with the iterative update formula to determine whether the objective function converges. If the objective function converges, the result is output. If the objective function does not converge, the next iteration is performed until the objective function converges.

2. The full waveform inversion method based on the self-adjoint VTI equation according to claim 1 is characterized in that: The specific formula of the split VTI acoustic wave equation is as follows: Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represents the normal stress in the horizontal and vertical directions respectively; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent the density, P-wave velocity and anisotropy parameter, respectively.

3. The full waveform inversion method based on the self-adjoint VTI equation according to claim 2 is characterized in that: The specific formula of the split VTI elastic wave equation is as follows: Among them, τ xx =τ xx1 +τ xx2 , τ xx and τ zz Represent the normal stress in the horizontal and vertical directions respectively; τ xz represents the shear stress in the shear direction; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter, respectively.

4. The full waveform inversion method based on the self-adjoint VTI equation according to claim 3 is characterized in that: The specific formula of the self-adjoint VTI acoustic wave equation is as follows: Among them, p and q represent auxiliary variables; τ xx2 represents the normal stress in the horizontal direction after splitting; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , ε and δ represent density, P-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

5. The full waveform inversion method based on the self-adjoint VTI equation according to claim 4 is characterized in that: The specific formula of the self-adjoint VTI elastic wave equation is as follows: Among them, p and q represent auxiliary variables; τ xx2 represents the normal stress in the horizontal direction after splitting; v x and v z Respectively represent the particle vibration speed in the horizontal and vertical directions; ρ, v p , v s , ε and δ represent density, P-wave velocity, S-wave velocity and anisotropy parameter respectively; s represents the earthquake source; and They represent the spatial first-order derivative operator in the horizontal direction, the spatial first-order derivative operator in the vertical direction, and the time first-order derivative operator respectively.

6. The full waveform inversion method based on the self-adjoint VTI equation according to claim 5 is characterized in that: The gradient formula of the objective function to the model parameters is as follows: Where "T" represents the transpose of the matrix; m represents the model parameter; u represents the forward wave field, A represents the forward operator; u * represents the accompanying wave field; x represents the spatial coordinate of a point underground.

7. The full waveform inversion method based on the self-adjoint VTI equation according to claim 6 is characterized in that: The objective function formula is as follows: Among them, m represents the model parameter, d pre represents the predicted data generated by model m, d obs Represents observation data collected in the field.

8. A full waveform inversion device based on a self-adjoint VTI equation, characterized in that: include: The first processing unit is used to process the normal stress τ in the VTI acoustic wave equation and the VTI elastic wave equation. xx Split into τ xx1 and τ xx2 , and obtain the split VTI acoustic wave equation and VTI elastic wave equation; The second processing unit is used to construct auxiliary variables p and q and bring them into the split VTI acoustic wave equation and VTI elastic wave equation to obtain a self-accompanied VTI acoustic wave equation and a self-accompanied VTI elastic wave equation; The third processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and calculate the forward wave field u and the predicted data d pre ; The fourth processing unit is used to use the objective function and the prediction data d pre Calculate the residual d res ; The fifth processing unit is used to solve the self-accompanied VTI acoustic wave equation and the self-accompanied VTI elastic wave equation, and back propagate the residual d res Calculate the accompanying wave field u * ; The sixth processing unit is used to use the forward wave field u and the accompanying wave field u based on the gradient formula of the objective function to the model parameters * Computing Gradients The seventh processing unit is used to calculate the step size α and pass the gradient The model m is updated with the iterative update formula to determine whether the objective function converges. If the objective function converges, the result is output. If the objective function does not converge, the next iteration is performed until the objective function converges.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the steps of the full waveform inversion method based on the self-adjoint VTI equation described in any one of claims 1 to 7 are implemented.

10. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the steps of the full waveform inversion method based on the self-adjoint VTI equation described in any one of claims 1 to 7 are implemented.

Images