Seismic wave impedance inversion method and system based on minimum and maximum concave high-order total variation constraint
By using the minimum-maximum concavity high-order total variation constraint method, the problem of insufficient inversion accuracy in the existing technology is solved, and higher inversion accuracy and resolution are achieved, artifacts are reduced, and the accuracy and continuity of the inversion results are improved.
Patent Information
- Application Number
- CN202410934630.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2024-07-12
- Publication Date
- 2025-01-10
- Estimated Expiration
- 2044-07-12
AI Technical Summary
Existing seismic impedance inversion methods based on total variation constraints have low inversion accuracy in structurally complex regions and produce artifacts in boundary regions, resulting in insufficient inversion accuracy.
The minimum-maximum concavity higher-order total variation constraint method is adopted. By constructing a linear relationship between the reflection coefficient and the natural logarithm of the wave impedance, and combining the minimum-maximum concavity higher-order total variation constraint term and the initial model constraint term, the logarithm of the wave impedance is iteratively updated using the alternating direction multiplier method. The augmented Lagrangian function is constructed and transformed into an unconstrained optimization problem.
It improves the accuracy and resolution of the inversion results, reduces artifacts in the boundary region, and enhances the continuity and accuracy of the inversion results. In particular, it significantly improves the clarity and reliability of the inversion results under complex geological conditions.
Smart Images

Figure CN119291771A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical field of seismic inversion and oil and gas reservoir prediction, and in particular to a seismic wave impedance inversion method and system using minimum-maximum-concave high-order total variation constraints. Background Art
[0002] Seismic impedance inversion is an important method for predicting oil and gas reservoirs. It establishes forward equations and objective functions through the mathematical and physical relationship between post-stack seismic records and seismic impedance, and uses optimization methods to solve the objective function to obtain the optimal solution of seismic impedance. Seismic impedance inversion technology based on sparse constraints is an important method of seismic inversion. It can improve the accuracy of inversion results by introducing sparse constraints into the forward equation.
[0003] The seismic impedance inversion method based on the traditional total variation constraint is an important method of the seismic impedance inversion technology based on sparse constraints. This method uses the first-order total variation constraint term to construct the forward model and objective function, and iterates the solution through the alternating direction multiplier method to finally obtain the impedance inversion result. This method can improve the inversion accuracy. However, this method has a weak ability to describe sparsity, resulting in low inversion accuracy in complex structural areas and artifacts in boundary areas, thereby reducing the inversion accuracy. Summary of the invention
[0004] In view of the defects of the prior art, the present invention provides a minimum-maximum concave high-order total variation constrained seismic wave impedance inversion method and system.
[0005] In order to achieve the above invention object, the technical solution adopted by the present invention is as follows:
[0006] A minimum-maximum concave high-order total variation constrained seismic wave impedance inversion method comprises the following steps:
[0007] S1. Preprocess the seismic data to obtain initial data for seismic inversion, including: limited angle or post-stack seismic data, impedance data corresponding to the target formation segment, seismic wavelet data of the target layer, and initial seismic inversion model;
[0008] S2. Construct the objective function of seismic wave impedance inversion constrained by minimum-maximum concave high-order total variation constraints, including: constructing a linear relationship model between the reflection coefficient and the natural logarithm of the wave impedance based on the data fidelity term, constructing a sparse constraint term based on the minimum-maximum concave high-order total variation constraint term, and constructing the objective function of seismic wave impedance inversion based on the initial model constraint term.
[0009] S3. Use the alternating direction multiplier method to iteratively update the wave impedance logarithm, calculate and update the natural logarithm of the wave impedance, and update the Lagrange multiplier terms and dual terms.
[0010] S4. Calculate the inversion result of seismic wave impedance according to the conditions, and determine whether the value of the natural logarithm of the wave impedance before and after the update meets the iterative update conditions. If so, use the updated natural logarithm of the wave impedance, the Lagrange multiplier term and the dual term as the initial values and transfer to the inversion iterative update module. If not, perform exponential operation on the updated natural logarithm of the wave impedance to obtain the inverted wave impedance result.
[0011] S5. Output and store the final inversion wave impedance results for subsequent analysis and application.
[0012] Furthermore, the step S1 includes the following sub-steps:
[0013] S101. Perform stacking processing on the pre-stack seismic data according to the incident angle to obtain limited-angle seismic data or post-stack seismic data, referred to as seismic data, which is denoted by S;
[0014] S102. Obtain impedance data corresponding to the target formation segment using logging data, where the impedance at the single-channel logging point is recorded as Z0;
[0015] S103. Using the well side channel seismic data and the target layer impedance at the logging location, obtain the target layer seismic wavelet data w, w = [w1, w2, ..., w p ], where the length of the seismic wavelet data is p, w i (i=1,2,……,p) is the wavelet data value at the i-th sampling moment;
[0016] S104. According to the target stratum layer data, the well logging data is low-pass filtered and interpolated according to the layer to obtain the initial model section of the wave impedance of the target layer segment.
[0017] Furthermore, the step S2 includes the following sub-steps:
[0018] S201. Calculate the initial natural logarithm of the wave impedance L0. The natural logarithm calculation formula is:
[0019] L0=ln(Z0)
[0020] Among them, Z0 represents the initial wave impedance, L0 represents the natural logarithmic wave impedance, and ln represents the natural logarithmic operation formula;
[0021] S202. Construct a linear calculation model of the reflection coefficient R and its wave impedance natural logarithm L, the calculation formula is:
[0022] R=DL
[0023] In the above formula, R is the reflection coefficient, L is the natural logarithm of the wave impedance, and D is the difference matrix, where
[0024] L=ln(Z)
[0025]
[0026] m represents the total sampling points of the wave impedance logarithm L;
[0027] S203. Construct the minimum maximum concave high-order total variation constraint. The minimum maximum concave high-order total variation constraint is calculated as:
[0028]
[0029] Among them, D xx represents the transverse second-order difference, D yy represents the longitudinal second-order difference, γ represents the minimum and maximum concave smoothing constraint parameter, whose value needs to be estimated in advance, ν represents the Morrow envelope function factor, |||·|||2 represents the L2 norm calculation formula, |||·|||1 represents the L1 norm calculation formula, where D xx L and LD yy The calculation formulas are:
[0030] (D xx L) s,t =2L s,t -L s-1,t -L s+1,t
[0031] (LD yy ) s,t =2L s,t -L s,t-1 -L s,t+1
[0032] S204. Construct a seismic impedance forward model using minimum-maximum concave high-order total variation constraints:
[0033]
[0034] In the above formula, S represents seismic data, W represents the seismic wavelet convolution matrix, which is composed of the seismic wavelet data w in step S103, the calculation method of the minimum and maximum concave high-order total variation constraint term is the same as that in step S203, λ represents the sparse constraint weight coefficient, whose value needs to be specified in advance, α represents the initial model constraint weight coefficient, whose value needs to be specified in advance, wherein the structure of the wavelet convolution matrix W is:
[0035]
[0036] p represents the length of the seismic wavelet w.
[0037] S205. Construct an augmented Lagrangian function in the seismic wave impedance forward model to transform the constrained optimization problem into an unconstrained optimization problem:
[0038]
[0039] In the above formula, R xx and R yy represents the Lagrange multiplier, C xx and C yy is the dual term, and β is a parameter greater than 0.
[0040] Furthermore, the step S3 includes the following sub-steps:
[0041] S301. Using the initial model wave impedance logarithm L0 as the wave impedance natural logarithm to iteratively calculate the initial value L 0 , take the initial value of Lagrange multiplier R xx =R yy = 0 and the initial value of the dual term C xx =C yy =0;
[0042] S302. Update the natural logarithm of the wave impedance L using the solution of the Sylvester equation k+1 , the updated calculation formula is:
[0043]
[0044] Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, vec represents the conversion of the matrix into a column vector, E represents the unit matrix, and T represents the transpose of the matrix. represents the Kronecker product operation, where A, B, and C are expressed as follows:
[0045]
[0046] S303. Update the Lagrange multiplier term R xx and R yy , the updated calculation formula is:
[0047]
[0048] in, and represents the Lagrange multiplier result of the k+1th iteration, sign represents the sign function, and The calculation formula is:
[0049]
[0050] Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, and represents the Lagrange multiplier result of the kth iteration, and represents the dual term result of the kth iteration;
[0051] S304. Update the dual term C xx and C yy , the updated calculation formula is:
[0052]
[0053] in and represents the dual term result of the k+1th iteration, and represents the dual term result of the kth iteration, L k+1 represents the k+1th iteration of the wave impedance logarithm, and Represents the Lagrange multiplier result of the k+1th iteration.
[0054] Furthermore, the step S4 includes the following sub-steps:
[0055] S401. Update the logarithm of seismic wave impedance L of the k+1th iteration k+1 ;
[0056] S402. If || L k+1 -L k ||2 / ||L k ||2>ε, where ε=5×10 -5 , then return to step S302 to loop; otherwise, calculate the wave impedance Z, and the wave impedance Z calculation formula is:
[0057] Z = exp(L k+1 )
[0058] Here, exp represents the natural exponential operation.
[0059] The present invention also discloses a seismic wave impedance inversion system, which can be used to implement the above-mentioned seismic wave impedance inversion method, specifically comprising:
[0060] Data preprocessing module: used to preprocess seismic data, including obtaining post-stack seismic record data, well logging data, seismic wavelets and seismic inversion initial model.
[0061] Inversion objective function construction module: Create the seismic wave impedance inversion calculation objective function based on the data fidelity term, the minimum and maximum concave high-order total variation constraint term, and the initial model constraint term.
[0062] Inversion iterative update module: the natural logarithm of the wave impedance is updated according to the alternating direction multiplier method and the seismic wave impedance inversion calculation objective function to obtain the updated natural logarithm of the wave impedance.
[0063] Inversion result judgment module: judge whether the value of the natural logarithm of the wave impedance before and after the update meets the iterative update conditions. If so, the updated natural logarithm of the wave impedance, the Lagrange multiplier term and the dual term are used as the initial values and transferred to the inversion iterative update module. If not, the updated natural logarithm of the wave impedance is subjected to exponential operation to obtain the inverted wave impedance result.
[0064] Result output module: output and store the final inversion wave impedance results, including generating visual charts and detailed reports, ensuring that all output data and charts meet the expected format and accuracy requirements.
[0065] The present invention also discloses a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor implements the above-mentioned seismic wave impedance inversion method when executing the program.
[0066] The present invention also discloses a computer-readable storage medium on which a computer program is stored. When the program is executed by a processor, the above-mentioned seismic wave impedance inversion method is implemented.
[0067] Compared with the prior art, the advantages of the present invention are:
[0068] 1. The minimum-maximum-concave high-order total variation constraint method is adopted to describe high-frequency sparse information through the minimum-maximum-concave method, which can more accurately capture the subtle features and changes in the seismic data and improve the inversion accuracy.
[0069] 2. Using high-order total variation to describe lateral and vertical differential information can better maintain and reflect the spatial variation law of seismic data and enhance the continuity and accuracy of inversion results.
[0070] 3. By combining minimum-maximum concavity and high-order total variation constraints, the wave impedance can be reconstructed more accurately in areas with complex structures, significantly improving the resolution and accuracy of the inversion results.
[0071] 4. Since the minimum-maximum concave high-order total variation constraint can effectively suppress the step effect and reduce artifacts in the boundary area, the inversion result is more realistic and reliable.
[0072] 5. The present invention significantly reduces the error in the inversion results and improves the accuracy of the seismic wave impedance inversion results.
[0073] 6. In areas with complex impedance structures, the present invention can significantly improve the resolution, making the inversion results clearer and more accurate under complex geological conditions. BRIEF DESCRIPTION OF THE DRAWINGS
[0074] Figure 1 A flow chart of a seismic wave impedance inversion method according to an embodiment of the present invention;
[0075] Figure 2 is a signal processing flow chart of an embodiment of the present invention;
[0076] Figure 3 The impedance inversion results and error profiles of the marmousi2 model using the traditional total variation algorithm in the embodiment of the present invention, where (a) is the impedance profile and (b) is the error profile;
[0077] Figure 4 The impedance inversion results and error profiles of the marmousi2 model of the algorithm proposed in the embodiment of the present invention, where (a) is the impedance profile and (b) is the error profile. DETAILED DESCRIPTION
[0078] In order to make the purpose, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below with reference to the accompanying drawings and examples.
[0079] like Figure 1 As shown,
[0080] A seismic wave impedance inversion method using minimum-maximum concave high-order total variation constraints comprises the following steps:
[0081] S1. Preprocess the seismic data to obtain initial data for seismic inversion;
[0082] S101. Perform stacking processing on the pre-stack seismic data according to the incident angle to obtain limited-angle seismic data or post-stack seismic data (hereinafter referred to as seismic data), where the seismic data is denoted as S;
[0083] S102. Obtain impedance data corresponding to the target formation segment using logging data, where the impedance at the single-channel logging point is recorded as Z0;
[0084] S103. Using the well side channel seismic data and the target layer impedance at the logging location, obtain the target layer seismic wavelet data w, w = [w1, w2, ..., w p ], where the length of the seismic wavelet data is p, w i (i=1,2,……,p) is the wavelet data value at the i-th sampling moment;
[0085] S104. According to the target stratum data, the well logging data is extrapolated and interpolated according to the stratum and low-pass filtered to obtain the initial model section of the target layer impedance.
[0086] S2. Construct the seismic impedance inversion objective function constrained by the minimum-maximum concave high-order total variation constraint;
[0087] S201. Calculate the initial natural logarithm of the wave impedance L0. The natural logarithm calculation formula is:
[0088] L0=ln(Z0)
[0089] Among them, Z0 represents the initial wave impedance, L0 represents the natural logarithmic wave impedance, and ln represents the natural logarithmic operation formula;
[0090] S202. Construct a linear calculation model of the reflection coefficient R and its wave impedance natural logarithm L, the calculation formula is:
[0091] R=DL
[0092] In the above formula, R is the reflection coefficient, L is the natural logarithm of the wave impedance, and D is the difference matrix, where
[0093] L=ln(Z)
[0094]
[0095] m represents the total sampling points of the wave impedance logarithm L;
[0096] S203. Construct the minimum maximum concave high-order total variation constraint. The minimum maximum concave high-order total variation constraint is calculated as:
[0097]
[0098] Among them, D xx represents the transverse second-order difference, D yy represents the longitudinal second-order difference, γ represents the minimum and maximum concave smoothing constraint parameter, whose value needs to be estimated in advance, ν represents the Morrow envelope function factor, |||·|||2 represents the L2 norm calculation formula, |||·|||1 represents the L1 norm calculation formula, where D xx L and LD yy The calculation formulas are:
[0099] (D xx L) s,t =2L s,t -L s-1,t -L s+1,t
[0100] (LD yy ) s,t =2L s,t -L s,t-1 -L s,t+1
[0101] S204. Construct a seismic impedance forward model using minimum-maximum concave high-order total variation constraints:
[0102]
[0103] In the above formula, S represents seismic data, W represents the seismic wavelet convolution matrix, which is composed of the seismic wavelet data w in step S103, the calculation method of the minimum and maximum concave high-order total variation constraint term is the same as that in step S203, λ represents the sparse constraint weight coefficient, whose value needs to be specified in advance, α represents the initial model constraint weight coefficient, whose value needs to be specified in advance, wherein the structure of the wavelet convolution matrix W is:
[0104]
[0105] p represents the length of the seismic wavelet w.
[0106] S205. Construct an augmented Lagrangian function in the seismic wave impedance forward model to transform the constrained optimization problem into an unconstrained optimization problem:
[0107]
[0108] In the above formula, R xx and R yy represents the Lagrange multiplier, C xx and C yy is the dual term, and β is a parameter greater than 0.
[0109] S3. Iteratively update the wave impedance logarithm using the alternating direction multiplication method;
[0110] S301. Using the initial model wave impedance logarithm L0 as the wave impedance natural logarithm to iteratively calculate the initial value L 0 , take the initial value of Lagrange multiplier R xx =R yy = 0 and the initial value of the dual term C xx =C yy =0;
[0111] S302. Update the natural logarithm of the wave impedance L using the solution of the Sylvester equation k+1 , the updated calculation formula is:
[0112]
[0113] Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, vec represents the conversion of the matrix into a column vector, E represents the unit matrix, and T represents the transpose of the matrix. represents the Kronecker product operation, where A, B, and C are expressed as follows:
[0114]
[0115] S303. Update the Lagrange multiplier term R xx and R yy , the updated calculation formula is:
[0116]
[0117] in, and represents the Lagrange multiplier result of the k+1th iteration, sign represents the sign function, and The calculation formula is:
[0118]
[0119] Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, and represents the Lagrange multiplier result of the kth iteration, and represents the dual term result of the kth iteration;
[0120] S304. Update the dual term C xx and C yy , the updated calculation formula is:
[0121]
[0122] in and represents the dual term result of the k+1th iteration, and represents the dual term result of the kth iteration, L k+1 represents the k+1th iteration of the wave impedance logarithm, and Represents the Lagrange multiplier result of the k+1th iteration.
[0123] S4. Calculate the inversion results of seismic wave impedance according to the conditions.
[0124] S401. Update the logarithm of seismic wave impedance L of the k+1th iteration k+1 ;
[0125] S402. If || L k+1 -L k ||2 / ||L k ||2>ε, where ε=5×10 -5 , then return to step S302 to loop; otherwise, calculate the wave impedance Z, and the wave impedance Z calculation formula is:
[0126] Z = exp(L k+1 )
[0127] Here, exp represents the natural exponential operation.
[0128] like Figure 2 As shown, a seismic wave impedance inversion signal processing flow using minimum-maximum concave high-order total variation constraints includes the following steps:
[0129] Obtain the initial value of the wave impedance and calculate its natural logarithm;
[0130] Construct a wave impedance forward model with minimum-maximum concave high-order total variation constraints, which includes data fidelity terms, minimum-maximum concave high-order total variation constraints and initial model constraints;
[0131] The logarithm of the initial model wave impedance is used as the initial value of the iterative calculation of the natural logarithm of the wave impedance, and the initial values of the Lagrange multiplier term and the dual term are taken as 0;
[0132] Use the alternating direction multiplier method to update the wave impedance logarithm, Lagrange multiplier terms and dual terms;
[0133] It is determined whether the updated wave impedance logarithm meets the loop condition. If so, the updated wave impedance logarithm is used as the initial value; otherwise, the wave impedance is calculated according to the wave impedance logarithm.
[0134] In order to verify the effectiveness of the method proposed in the present invention, the method of the present invention is applied to the impedance inversion calculation of the marmousi2 model and compared with the inversion results of the traditional total variation (ATV) method. The inversion results are as follows: Figure 3 and Figure 4 As shown, Figure 3 (a) and Figure 3 (b) is the impedance inversion result profile and impedance error profile of the marmousi2 model based on the traditional total variation (ATV) method. Figure 4 (a) and Figure 4 (b) is the impedance inversion result profile and impedance error profile of the marmousi2 model of the proposed minimum-maximum-concave high-order total variation constraint algorithm. It can be seen from the figure that the inversion result error of the proposed method is significantly lower than that of the traditional method, and the resolution is significantly improved in areas with complex impedance structures, while suppressing artifacts in the boundary area.
[0135] In yet another embodiment of the present invention, a seismic wave impedance inversion system is provided. The system can be used to implement the above-mentioned seismic wave impedance inversion method, specifically comprising:
[0136] Data preprocessing module: used to preprocess seismic data, including obtaining post-stack seismic record data, well logging data, seismic wavelets and seismic inversion initial model.
[0137] Inversion objective function construction module: Create the seismic wave impedance inversion calculation objective function based on the data fidelity term, the minimum and maximum concave high-order total variation constraint term, and the initial model constraint term.
[0138] Inversion iterative update module: the natural logarithm of the wave impedance is updated according to the alternating direction multiplier method and the seismic wave impedance inversion calculation objective function to obtain the updated natural logarithm of the wave impedance.
[0139] Inversion result judgment module: judge whether the value of the natural logarithm of the wave impedance before and after the update meets the iterative update conditions. If so, the updated natural logarithm of the wave impedance, the Lagrange multiplier term and the dual term are used as the initial values and transferred to the inversion iterative update module. If not, the updated natural logarithm of the wave impedance is subjected to exponential operation to obtain the inverted wave impedance result.
[0140] Result output module: output and store the final inversion wave impedance results, including generating visual charts and detailed reports, ensuring that all output data and charts meet the expected format and accuracy requirements.
[0141] In another embodiment of the present invention, a terminal device is provided, the terminal device includes a processor and a memory, the memory is used to store a computer program, the computer program includes program instructions, and the processor is used to execute the program instructions stored in the computer storage medium. The processor may be a central processing unit (CPU), or other general-purpose processors, digital signal processors (DSP), application-specific integrated circuits (ASIC), field-programmable gate arrays (FPGA) or other programmable logic devices, discrete gates or transistor logic devices, discrete hardware components, etc., which are the computing core and control core of the terminal, which are suitable for implementing one or more instructions, and are specifically suitable for loading and executing one or more instructions to implement the corresponding method flow or corresponding functions; the processor described in the embodiment of the present invention can be used for the operation of the seismic wave impedance inversion method.
[0142] In another embodiment of the present invention, the present invention further provides a storage medium, specifically a computer-readable storage medium (Memory), which is a memory device in a terminal device for storing programs and data. It can be understood that the computer-readable storage medium here can include both the built-in storage medium in the terminal device and the extended storage medium supported by the terminal device. The computer-readable storage medium provides a storage space, which stores the operating system of the terminal. In addition, one or more instructions suitable for being loaded and executed by the processor are also stored in the storage space, and these instructions can be one or more computer programs (including program codes). It should be noted that the computer-readable storage medium here can be a high-speed RAM memory, or a non-volatile memory (non-volatile memory), such as at least one disk storage.
[0143] One or more instructions stored in a computer-readable storage medium can be loaded and executed by a processor to implement the corresponding steps of the seismic wave impedance inversion method in the above embodiment; one or more instructions in the computer-readable storage medium are loaded and executed by the processor.
[0144] Those skilled in the art will appreciate that embodiments of the present invention may be provided as methods, systems, or computer program products. Therefore, the present invention may take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.
[0145] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. 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 produce 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.
[0146] 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.
[0147] 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. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.
[0148] Those skilled in the art will appreciate that the embodiments described herein are intended to help readers understand the implementation methods of the present invention, and should be understood that the protection scope of the present invention is not limited to such special statements and embodiments. Those skilled in the art can make various other specific variations and combinations that do not deviate from the essence of the present invention based on the technical revelations disclosed by the present invention, and these variations and combinations are still within the protection scope of the present invention.
Claims
1. A minimum-maximum concave high-order total variation constrained seismic impedance inversion method, characterized in that: The following steps are involved: S1. Preprocess the seismic data to obtain initial data for seismic inversion, including: limited angle or post-stack seismic data, impedance data corresponding to the target formation segment, seismic wavelet data of the target layer, and initial seismic inversion model; S2. Constructing a seismic wave impedance inversion objective function constrained by a minimum-maximum concave high-order total variation constraint, including: constructing a linear relationship model between the reflection coefficient and the natural logarithm of the wave impedance according to a data fidelity term, constructing a sparse constraint term according to a minimum-maximum concave high-order total variation constraint term, and constructing an objective function of seismic wave impedance inversion according to an initial model constraint term; S3. Use the alternating direction multiplier method to iteratively update the wave impedance logarithm, calculate and update the natural logarithm of the wave impedance, and update the Lagrange multiplier term and the dual term; S4. Calculate the inversion result of seismic wave impedance according to the conditions, and judge whether the value of the natural logarithm of the wave impedance before and after the update meets the iterative update condition. If so, take the updated natural logarithm of the wave impedance, the Lagrange multiplier term and the dual term as the initial value and transfer to the inversion iterative update module. If not, perform exponential operation on the updated natural logarithm of the wave impedance to obtain the inversion wave impedance result; S5. Output and store the final inversion wave impedance results for subsequent analysis and application.
2. The seismic wave impedance inversion method according to claim 1, characterized in that: The step S1 comprises the following sub-steps: S101. Perform stacking processing on the pre-stack seismic data according to the incident angle to obtain limited-angle seismic data or post-stack seismic data, referred to as seismic data, which is denoted by S; S102. Obtain impedance data corresponding to the target formation segment using logging data, where the impedance at the single-channel logging point is recorded as Z0; S103. Using the well side channel seismic data and the target layer impedance at the logging location, obtain the target layer seismic wavelet data w, w = [w1, w2, ..., w p ], where the length of the seismic wavelet data is p, w i (i=1,2,……,p) is the wavelet data value at the i-th sampling moment; S104. According to the target stratum layer data, the well logging data is low-pass filtered and interpolated according to the layer to obtain the initial model section of the wave impedance of the target layer segment.
3. The seismic wave impedance inversion method according to claim 1, characterized in that: The step S2 comprises the following sub-steps: S201. Calculate the initial natural logarithm of the wave impedance L0. The natural logarithm calculation formula is: L0=ln(Z0) Among them, Z0 represents the initial wave impedance, L0 represents the natural logarithmic wave impedance, and ln represents the natural logarithmic operation formula; S202. Construct a linear calculation model of the reflection coefficient R and its wave impedance natural logarithm L, the calculation formula is: R=DL In the above formula, R is the reflection coefficient, L is the natural logarithm of the wave impedance, and D is the difference matrix, where L=ln(Z) m represents the total sampling points of the wave impedance logarithm L; S203. Construct the minimum maximum concave high-order total variation constraint. The minimum maximum concave high-order total variation constraint is calculated as: Among them, D xx represents the transverse second-order difference, D yy represents the longitudinal second-order difference, γ represents the minimum and maximum concave smoothing constraint parameter, whose value needs to be estimated in advance, ν represents the Morrow envelope function factor, |||·|||2 represents the L2 norm calculation formula, |||·|||1 represents the L1 norm calculation formula, where D xx L and LD yy The calculation formulas are: (D xx L) s,t =2L s,t -L s-1,t -L s+1,t (LD yy ) s,t =2L s,t -L s,t-1 -L s,t+1 S204. Construct a seismic impedance forward model using minimum-maximum concave high-order total variation constraints: In the above formula, S represents seismic data, W represents the seismic wavelet convolution matrix, which is composed of the seismic wavelet data w in step S103, the calculation method of the minimum and maximum concave high-order total variation constraint term is the same as that in step S203, λ represents the sparse constraint weight coefficient, whose value needs to be specified in advance, α represents the initial model constraint weight coefficient, whose value needs to be specified in advance, wherein the structure of the wavelet convolution matrix W is: p represents the length of the seismic wavelet w; S205. Construct an augmented Lagrangian function in the seismic wave impedance forward model to transform the constrained optimization problem into an unconstrained optimization problem: In the above formula, R xx and R yy represents the Lagrange multiplier, C xx and C yy is the dual term, and β is a parameter greater than 0.
4. The seismic wave impedance inversion method according to claim 1, characterized in that: The step S3 comprises the following sub-steps: S301. Using the initial model wave impedance logarithm L0 as the wave impedance natural logarithm to iteratively calculate the initial value L 0 , take the initial value of Lagrange multiplier R xx =R yy = 0 and the initial value of the dual term C xx =C yy =0; S302. Update the natural logarithm of the wave impedance L using the solution of the Sylvester equation k+1 , the updated calculation formula is: Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, vec represents the conversion of the matrix into a column vector, E represents the unit matrix, and T represents the transpose of the matrix. represents the Kronecker product operation, where A, B, and C are expressed as follows: S303. Update the Lagrange multiplier term R xx and R yy , the updated calculation formula is: in, and represents the Lagrange multiplier result of the k+1th iteration, sign represents the sign function, and The calculation formula is: Among them, L k+1 represents the k+1th iteration of the wave impedance logarithm, and represents the Lagrange multiplier result of the kth iteration, and represents the dual term result of the kth iteration; S304. Update the dual term C xx and C yy , the updated calculation formula is: in and represents the dual term result of the k+1th iteration, and represents the dual term result of the kth iteration, L k+1 represents the k+1th iteration of the wave impedance logarithm, and Represents the Lagrange multiplier result of the k+1th iteration.
5. The seismic wave impedance inversion method according to claim 1, characterized in that: The step S4 comprises the following sub-steps: S401. Update the logarithm of seismic wave impedance L of the k+1th iteration k+1 ; S402. If || L k+1 -L k ||2 / ||L k ||2>ε, where ε=5×10 -5 , then return to step S302 to loop; otherwise, calculate the wave impedance Z, and the wave impedance Z calculation formula is: Z=exp(L k+1 ) Here, exp represents the natural exponential operation.
6. A seismic wave impedance inversion system, characterized in that: The system can be used to implement the seismic wave impedance inversion method described in any one of claims 1 to 5, specifically comprising: Data preprocessing module: used to preprocess seismic data, including obtaining post-stack seismic record data, well logging data, seismic wavelets and seismic inversion initial model; Inversion objective function construction module: Create seismic wave impedance inversion calculation objective function based on data fidelity terms, minimum and maximum concave high-order total variation constraints, and initial model constraints; Inversion iterative update module: updates the natural logarithm of wave impedance according to the alternating direction multiplier method and the seismic wave impedance inversion calculation objective function to obtain the updated natural logarithm of wave impedance; Inversion result judgment module: judge whether the value of the natural logarithm of the wave impedance before and after the update meets the iterative update conditions. If so, the updated natural logarithm of the wave impedance, the Lagrange multiplier term and the dual term are used as the initial values and transferred to the inversion iterative update module. If not, the updated natural logarithm of the wave impedance is subjected to exponential operation to obtain the inverted wave impedance result. Result output module: output and store the final inversion wave impedance results, including generating visual charts and detailed reports, ensuring that all output data and charts meet the expected format and accuracy requirements.
7. A computer device, characterized in that: The invention comprises a memory, a processor and a computer program stored in the memory and executable on the processor, wherein when the processor executes the program, the seismic wave impedance inversion method as claimed in any one of claims 1 to 5 is implemented.
8. A computer-readable storage medium, characterized in that: A computer program is stored thereon, and when the program is executed by a processor, the seismic wave impedance inversion method described in one of claims 1 to 5 is implemented.
Citation Information
Patent Citations
Three-dimensional two-way acoustic wave equation pre-stack imaging systems and methods
CN101017204A
Seismic inversion method and system based on generalized total variation regularization
CN108037531A
Seismic wave impedance inversion method based on reweighted L1 norm sparse constraint
CN113640871A
Seismic wave impedance inversion method based on sparse constraint of non-convex arc tangent function zeta
CN114994757A
Seismic wave impedance inversion method and system based on logarithm total variation sparse constraint
CN115494547A
Cited By
Wave impedance inversion method based on reweighted arc tangent total variation constraint
CN119882044A
Phased post-stack seismic inversion method and system based on geologic body unit parameter updating
CN120143256A