CN114994757B - Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint - Google Patents

Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint Download PDF

Info

Publication number
CN114994757B
CN114994757B CN202210719263.2A CN202210719263A CN114994757B CN 114994757 B CN114994757 B CN 114994757B CN 202210719263 A CN202210719263 A CN 202210719263A CN 114994757 B CN114994757 B CN 114994757B
Authority
CN
China
Prior art keywords
wave impedance
representing
log
formula
seismic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210719263.2A
Other languages
Chinese (zh)
Other versions
CN114994757A (en
Inventor
吴昊
何良盛
文晓涛
陈禹
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Chengdu Univeristy of Technology
Original Assignee
Chengdu Univeristy of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Chengdu Univeristy of Technology filed Critical Chengdu Univeristy of Technology
Priority to CN202210719263.2A priority Critical patent/CN114994757B/en
Publication of CN114994757A publication Critical patent/CN114994757A/en
Application granted granted Critical
Publication of CN114994757B publication Critical patent/CN114994757B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/30Analysis
    • G01V1/306Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/62Physical property of subsurface
    • G01V2210/622Velocity, density or impedance
    • G01V2210/6226Impedance

Abstract

Disclosure of the inventionA seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint is disclosed, which comprises the following steps: step 1: preprocessing the seismic data and acquiring an initial model of wave impedance; and 2, step: calculating an initial log wave impedance L from the processed seismic data 0 Constructing a non-convex arc tangent function zeta, and establishing a forward model and a target function based on sparse constraint of the non-convex arc tangent function zeta; and 3, step 3: updating the log-wave impedance L by using an alternative direction multiplier method and a forward model to obtain the updated log-wave impedance L i (ii) a And 4, step 4: judging whether the values before and after updating satisfy L i+1 ‑L i || 2 /||L i || 2 If yes, returning to the step 3 for circulation; if not, performing exponential operation to obtain an inversion result of the wave impedance. The method eliminates the pseudo-layer phenomenon caused by L1 norm sparse constraint, and achieves the effect of improving the accuracy of the inversion result.

Description

Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint
Technical Field
The invention relates to the field of oil-gas exploration and seismic inversion, in particular to a seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint.
Background
The inversion of seismic wave impedance obtains the process of seismic wave impedance from seismic records by an optimization method according to the mathematical physical relationship between the post-stack seismic records and the wave impedance, and is an important means of oil-gas exploration. The seismic wave impedance inversion technology based on sparse constraint is an important method for seismic wave impedance inversion. Sparse constraint is used as prior information and introduced into the forward model, and the accuracy of an inversion result is successfully improved.
The seismic wave impedance inversion method based on the L1 norm sparse constraint is an important method of the seismic wave impedance inversion technology based on the sparse constraint. The method adopts L1 norm as a constraint term to construct a forward model, and an optimization problem formed by the forward model is iteratively solved by an optimization method to finally obtain a wave impedance inversion result. Tayler et al used the L1 norm for seismic inversion for the first time and demonstrated that the method can improve inversion accuracy. Wang et al used the L1 norm constraint for single-pass wave impedance inversion and solved using a gradient descent method, demonstrating the feasibility of the method in seismic wave impedance inversion. However, the L1 norm cannot completely mine the sparsity, so that a pseudo layer appears in an inversion result, and the precision is reduced; there is therefore a need for an inversion method that overcomes the above problems.
Disclosure of Invention
The invention aims to: the invention provides a seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint, solves the problem that a pseudo layer caused by an L1 norm is adopted in the existing seismic wave impedance inversion method based on L1 norm sparse constraint, and achieves the effect of improving the accuracy of inversion results.
The technical scheme adopted by the invention is as follows:
a seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint comprises the following steps:
step 1: preprocessing the seismic data and acquiring an initial model of wave impedance;
step 2: calculating an initial log wave impedance L from the processed seismic data 0 Constructing a non-convex arc tangent function zeta, and establishing a forward model and a target function based on sparse constraint of the non-convex arc tangent function zeta;
and step 3: updating the log-wave impedance L by using an alternative direction multiplier method and a forward model to obtain the updated log-wave impedance L i
And 4, step 4: judging whether the values before and after updating satisfy L i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3 for circulation; if not, performing exponential operation to obtain an inversion result of the wave impedance;
wherein L is i+1 Represents the i +1 th iteration result of the log-wave impedance L i The result of the i-th iteration of the log-wave impedance L is represented and tol represents the error.
Preferably, the step 1 comprises the following steps:
step 1.1: input post-stack seismic record S 0 Wavelet data w and logging data;
step 1.2: and extracting the horizon information and carrying out interpolation filtering on the logging data to obtain an initial model of the wave impedance of the parameter to be inverted.
Preferably, the step 2 comprises the following steps:
step 2.1: obtaining initial logarithmic wave impedance L by an initial model based on wave impedance according to a logarithmic operation formula 0 (ii) a The logarithmic operation formula is shown in formula (1):
l = ln (Z) formula (1)
Wherein Z represents wave impedance, L represents logarithmic wave impedance, and ln represents logarithmic operation;
step 2.2: constructing a non-convex arc tangent function ζ as shown in equation (2):
Figure BDA0003709821340000021
wherein x represents data to be constrained, a represents the non-convex degree of an arctangent function zeta, | | represents an absolute value, and arctan represents arctangent;
step 2.3: constructing a seismic wave impedance forward modeling based on a non-convex arctangent function zeta sparse constraint, wherein the modeling is as shown in a formula (3), a formula (4) and a formula (5):
Figure BDA0003709821340000022
Figure BDA0003709821340000023
Figure BDA0003709821340000031
wherein λ represents a sparse constraint weight coefficient, β represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, | | | | calving F Represents Frobenius norm, m represents the total sampling point of log-wave impedance, L n Representing the value of the nth sample point, w, in the log-wave impedance L 1 1 st sample point, w, representing a seismic wavelet q Q is the number of the q sampling points of the seismic wavelet, q is the total number of sampling points of the seismic wavelet, L 0 Representing the initial log-wave impedance, a representing the degree of non-convexity of the arctangent function ζ, S 0 Representing a post-stack seismic record;
step 2.4: introducing a Lagrange multiplier term R and a dual term C on the basis of a forward model to obtain a seismic wave impedance inversion target function based on the non-convex arctangent function zeta sparse constraint, as shown in a formula (6):
Figure BDA0003709821340000032
wherein gamma is a dual term weight coefficient, gamma/lambda represents gamma divided by lambda, lambda represents a sparse constraint weight coefficient, beta represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, and | | Y F Represents Frobenius norm, m represents the total sampling point of log-wave impedance, L represents log-wave impedance, L 0 Representing the initial log-wave impedance, a representing the degree of non-convexity of the arctangent function ζ, S 0 Representing post-stack seismic records, R n Representing the nth value in the lagrange multiplier term R.
Preferably, the step 3 comprises the following steps:
step 3.1: impedance L of initial log wave 0 Initial model L as a sparsity constraint for a non-convex arctangent function ζ 1 And introducing an initial value R of a Lagrange multiplier term 1 =0 and its dual initial value C 1 =0;
Step 3.2: updating wave impedance logarithm L by using alternative direction multiplier algorithm and criterion of taking extreme value when gradient is zero i +1 The calculation formula is shown in formula (7):
Figure BDA0003709821340000033
wherein E represents an identity matrix, D T Representing the transpose of the difference matrix D, W T Representing the transpose, L, of the wavelet convolution matrix W i+1 The result of the (i + 1) th iteration, R, representing the log-wave impedance L i Representing the result of the i-th iteration of the Lagrange multiplier term R, C i Representing the ith iteration result of the even term C, T representing the transposition of the matrix, i representing the ith iteration, gamma being the weight coefficient of the even term, beta representing the weight coefficient of the initial model, L 0 Representing the initial log-wave impedance, S 0 Representing a post-stack seismic record;
step 3.3: updating Lagrange multiplier term R according to alternative direction multiplier algorithm and high-order approximation method i+1 The calculation formula is shown as formula (8):
Figure BDA0003709821340000041
wherein the symbol/represents a dot division operation,
Figure BDA0003709821340000042
representing a matrix dot product operation, a representing the degree of non-convexity of the arctangent function ζ, D representing a difference matrix, C i Represents the result of the i-th iteration on the even term C, L i+1 The result of the (i + 1) th iteration, R, representing the log-wave impedance L i+1 Representing the (i + 1) th iteration result of the Lagrange multiplier term R;
sign represents a sign function, and the expression is shown as formula (9):
Figure BDA0003709821340000043
wherein x is j The expression vector DL i+1 +C i The jth element of (1);
f 3 and f 2 For intermediate variables, the expressions are shown in equation (10) and equation (11):
f 3 =(sqrt(f 1 .^2-f 2 .^3)-f 1 ) A 1/3 formula (10)
f 2 =b.^2./(9.*a 2 )-b./(3.*a 2 ) Formula (11)
Wherein, the symbol ^ represents a point power operation, sqrt represents a square root function, and a represents the non-convexity degree of an arctangent function ζ;
f 1 and b is expressed as shown in equation (12) and equation (13):
Figure BDA0003709821340000044
b=1-a*abs(R i ) Formula (13)
Wherein R is i Representing the ith iteration result of the Lagrange multiplier term R, wherein a represents the non-convex degree of the arctangent function zeta; gamma is a dual term weight coefficient, lambda represents a sparse constraint weight coefficient, and abs represents an absolute value function;
step 3.4: updating the parity term C according to the alternative direction multiplier algorithm and Fermat theorem i+1 The calculation formula is shown as formula (14):
C i+1 =C i +(DL i+1 -R i+1 ) Formula (14)
Wherein, C i+1 Represents the i +1 th iteration result to the even term C, D represents the difference matrix, C i Represents the result of the i-th iteration on the even term C, L i+1 Represents the (i + 1) th iteration result of the log-wave impedance L; r i+1 Represents the i +1 th iteration result of the Lagrangian multiplier term R.
Preferably, the step 4 determines whether | | | L is satisfied i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3.2 for circulation, and if not, calculating the wave impedance Z according to exponential operation r As shown in equation (15):
Figure BDA0003709821340000051
wherein the content of the first and second substances,
Figure BDA0003709821340000052
denotes the natural logarithm e as the base number, L i+1 Is an index, L i+1 Represents the i +1 th iteration result of the log-wave impedance L.
In summary, due to the adoption of the technical scheme, the invention has the beneficial effects that:
according to the method, the L1 norm in the L1 norm sparse constraint is replaced by the non-convex arc tangent function zeta, the pseudo layer phenomenon caused by the L1 norm sparse constraint is successfully eliminated by introducing the non-convex arc tangent function zeta, and the effect of improving the accuracy of the inversion result is achieved.
Drawings
The invention will now be described, by way of example, with reference to the accompanying drawings, in which:
FIG. 1 is a block diagram of a flow chart of a seismic wave impedance inversion method based on a non-convex arctangent function zeta sparse constraint according to the invention;
FIG. 2 is a noisy seismic section of the present invention;
FIG. 3 is a diagram of an initial model of multi-channel impedance in accordance with the present invention;
FIG. 4 is a schematic diagram of the inversion result of multi-channel impedance of the present invention;
FIG. 5 is a diagram illustrating the updating result and comparison of the bypass wave impedance of the single well according to the present invention.
Detailed Description
In order to make the objects, technical solutions and advantages of the present invention more apparent, the present invention is further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the detailed description and specific examples, while indicating the preferred embodiment of the invention, are intended for purposes of illustration only and are not intended to limit the scope of the invention. The components of embodiments of the present invention generally described and illustrated in the figures herein may be arranged and designed in a wide variety of different configurations.
Thus, the following detailed description of the embodiments of the present invention, presented in the figures, is not intended to limit the scope of the invention, as claimed, but is merely representative of selected embodiments of the invention. All other embodiments, which can be derived by a person skilled in the art from the embodiments of the present invention without making any creative effort, shall fall within the protection scope of the present invention.
It is noted that relational terms such as "first" and "second," and the like, may be used solely to distinguish one entity or action from another entity or action without necessarily requiring or implying any actual such relationship or order between such entities or actions. Also, the terms "comprises," "comprising," or any other variation thereof, are intended to cover a non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process, method, article, or apparatus. Without further limitation, an element defined by the phrases "comprising one of 8230; \8230;" 8230; "does not exclude the presence of additional like elements in a process, method, article, or apparatus that comprises the element.
The present invention is described in detail below with reference to fig. 1-5.
Example 1
The technical problem to be solved is as follows: the pseudo-layer phenomenon caused by L1 norm sparse constraint is eliminated, and the effect of improving the accuracy of the inversion result is achieved.
The technical means is as follows: a seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint comprises the following steps:
step 1: preprocessing the seismic data and acquiring an initial model of wave impedance;
step 2: calculating an initial log wave impedance L from the processed seismic data 0 Constructing a non-convex arc tangent function zeta, and establishing a forward model and a target function based on sparse constraint of the non-convex arc tangent function zeta;
and 3, step 3: updating the log-wave impedance L by using an alternative direction multiplier method and a forward model to obtain the updated log-wave impedance L i
And 4, step 4: judging whether the values before and after updating satisfy L i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3 for circulation; if not, performing exponential operation to obtain an inversion result of the wave impedance;
wherein L is i+1 The i +1 th iteration result representing the log-wave impedance L, L i The result of the i-th iteration of the log-wave impedance L is represented and tol represents the error.
Example 2
The step 1 comprises the following steps:
step 1.1: input post-stack seismic record S 0 Wavelet data w and logging data;
step 1.2: and extracting the horizon information and carrying out interpolation filtering on the logging data to obtain an initial model of the wave impedance of the parameter to be inverted.
Example 3
The step 2 comprises the following steps:
step 2.1: obtaining initial logarithmic wave impedance L by an initial model based on wave impedance according to a logarithmic operation formula 0 (ii) a The logarithmic operation formula is shown in formula (1):
l = ln (Z) formula (1)
Wherein Z represents wave impedance, L represents logarithmic wave impedance, ln represents logarithmic operation;
step 2.2: constructing a non-convex arctangent function ζ as shown in equation (2):
Figure BDA0003709821340000071
wherein, x represents the data to be constrained, a represents the non-convex degree of an arctan function zeta, | represents an absolute value, and arctan represents arctan;
step 2.3: constructing a seismic wave impedance forward modeling based on a non-convex arctangent function zeta sparse constraint, wherein the modeling is as shown in a formula (3), a formula (4) and a formula (5):
Figure BDA0003709821340000072
Figure BDA0003709821340000073
Figure BDA0003709821340000074
wherein λ represents a sparse constraint weight coefficient, β represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, | | | | | survival F Representing the Frobenius norm, m representing the total sampling point of the log-wave impedance, L n Representing the value of the nth sample point, w, in the log-wave impedance L 1 1 st sample point, w, representing a seismic wavelet q Represents the q-th sampling point of the seismic wavelet, q represents the total sampling point number of the seismic wavelet, L 0 Denotes the initial log wave impedance, a denotes the degree of non-convexity of the arctangent function ζ, S 0 Representing a post-stack seismic record;
step 2.4: a Lagrange multiplier term R and a dual term C are introduced on the basis of a forward model, and a seismic wave impedance inversion target function based on a non-convex arctangent function zeta sparse constraint is obtained and is shown in a formula (6):
Figure BDA0003709821340000081
wherein gamma is a dual term weight coefficient, gamma/lambda represents gamma divided by lambda, lambda represents a sparse constraint weight coefficient, beta represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, and | | | | y F Representing the Frobenius norm, m representing the total sampling point of the log-wave impedance, L representing the log-wave impedance, L 0 Representing the initial log-wave impedance, a representing the degree of non-convexity of the arctangent function ζ, S 0 Representing post-stack seismic records, R n Representing the 1 nth value in the lagrange multiplier term R.
Example 4
The step 3 comprises the following steps:
step 3.1: impedance L of initial log wave 0 As non-convex reversalInitial model L of tangent function zeta sparse constraint 1 Introducing initial value R of Lagrange multiplier term 1 =0 and its dual initial value C 1 =0;
Step 3.2: updating wave impedance logarithm L by using alternative direction multiplier algorithm and criterion of zero gradient extremum i +1 The calculation principle is shown in formula (7):
L i+1 =(D T W T WD+βE+γD T D) -1 [D T W T S 0 +βL 0 +γD T (R i -C i )]formula (7)
Wherein E represents an identity matrix, D T Representing the transpose of a difference matrix D, W T Representing the transpose, L, of a wavelet convolution matrix W i+1 The result of the (i + 1) th iteration, R, representing the log-wave impedance L i Representing the result of the i-th iteration of the Lagrange multiplier term R, C i Representing the ith iteration result of the even term C, T representing the transposition of the matrix, i representing the ith iteration, gamma being the weight coefficient of the even term, beta representing the weight coefficient of the initial model, L 0 Representing the initial log-wave impedance, S 0 Representing a post-stack seismic record;
step 3.3: updating Lagrange multiplier term R according to alternative direction multiplier algorithm and high-order approximation method i+1 The calculation formula is shown in formula (8):
Figure BDA0003709821340000091
wherein the symbol/represents a dot division operation,
Figure BDA0003709821340000092
representing a matrix dot product operation, a representing the degree of non-convexity of the arctangent function ζ, D representing a difference matrix, C i Denotes the result of the i-th iteration, L, on the even term C i+1 Represents the i +1 th iteration result of the log-wave impedance L, R i+1 Representing the i +1 th iteration result of the Lagrange multiplier term R;
sign represents a sign function, and the expression is shown as formula (9):
Figure BDA0003709821340000093
wherein x is j The expression vector DL i+1 +C i The jth element of (1);
f 3 and f 2 For intermediate variables, the expressions are as shown in equations (10) and (11):
f 3 =(sqrt(f 1 .^2-f 2 .^3)-f 1 ) A 1/3 formula (10)
f 2 =b.^2./(9.*a 2 )-b./(3.*a 2 ) Formula (11)
Wherein, the symbol ^ represents a point power operation, sqrt represents a square root function, and a represents the non-convexity degree of an arctangent function ζ;
f 1 the expression of b is shown in equation (12) and equation (13):
Figure BDA0003709821340000094
b=1-a*abs(R i ) Formula (13)
Wherein R is i Representing the ith iteration result of the Lagrange multiplier term R, and a representing the non-convex degree of an arctangent function zeta; gamma is a dual term weight coefficient, lambda represents a sparse constraint weight coefficient, and abs represents an absolute value function;
step 3.4: updating the even term C according to the alternative direction multiplier algorithm and the Fermat lemma i+1 The calculation formula is shown as formula (14):
C i+1 =C i +(DL i+1 -R i+1 ) Formula (14)
Wherein, C i+1 Represents the i +1 th iteration result to the even term C, D represents the difference matrix, C i Denotes the result of the i-th iteration, L, on the even term C i+1 Represents the (i + 1) th iteration result of the log-wave impedance L; r is i+1 Representing LagrangeThe (i + 1) th iteration result of the daily multiplier term R.
Example 5
In step 4, it is determined whether | | | L is satisfied i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3.2 for circulation, and if not, calculating the wave impedance Z according to exponential operation r As shown in equation (15):
Figure BDA0003709821340000101
wherein the content of the first and second substances,
Figure BDA0003709821340000102
denotes the natural logarithm e as the base number, L i+1 Is an index, L i+1 Represents the i +1 th iteration result of the log-wave impedance L.
In conclusion, the L1 norm in the L1 norm sparse constraint is replaced by the non-convex arc tangent function zeta, the introduction of the non-convex arc tangent function zeta successfully eliminates the pseudo layer phenomenon caused by the L1 norm sparse constraint, and the effect of improving the accuracy of the inversion result is achieved.
And (3) effect analysis: 2-4, the inversion Result of the wave impedance profile can better follow the trend of seismic record compared with the Initial model, and the correctness of the method is proved, wherein an Initial model represents the wave impedance Initial model, result of arctan represents the updated wave impedance Result, true AI represents the True wave impedance, the inversion Result is very close to the True wave impedance, and no pseudo-layer phenomenon is caused, and the correctness of the method is further proved. The invention provides a seismic wave impedance inversion method based on sparse constraint of a non-convex arctangent function zeta by replacing an L1 norm in an L1 norm sparse constraint with the non-convex arctangent function zeta and combining an alternating direction multiplier method. And carrying out inversion updating on the initial model by adopting the inversion method, and outputting an optimal inversion result through repeated iteration. The seismic wave impedance inversion method based on the L1 norm sparse constraint solves the problem of a pseudo layer phenomenon caused by the adoption of the L1 norm in the existing seismic wave impedance inversion method based on the L1 norm sparse constraint, and achieves the effect of improving the accuracy of the inversion result.

Claims (5)

1. A seismic wave impedance inversion method based on non-convex arctangent function zeta sparse constraint is characterized by comprising the following steps:
step 1: preprocessing the seismic data and acquiring an initial model of wave impedance;
step 2: calculating an initial log wave impedance L from the processed seismic data 0 Constructing a non-convex arc tangent function zeta, and establishing a forward model and a target function based on sparse constraint of the non-convex arc tangent function zeta;
and 3, step 3: updating the log-wave impedance L by using an alternative direction multiplier method and a forward model to obtain the updated log-wave impedance L i
And 4, step 4: judging whether the values before and after updating satisfy L i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3 for circulation; if not, performing exponential operation to obtain an inversion result of the wave impedance;
wherein L is i+1 Represents the i +1 th iteration result of the log-wave impedance L i Represents the result of the ith iteration of the log wave impedance L, and tol represents the error.
2. The seismic wave impedance inversion method based on the sparsity constraint of the non-convex arctangent function ζ according to claim 1, wherein:
the step 1 comprises the following steps:
step 1.1: input post-stack seismic record S 0 Wavelet data w and well logging data;
step 1.2: and extracting horizon information and carrying out interpolation filtering on the logging data to obtain an initial model of the wave impedance of the parameter to be inverted.
3. The seismic wave impedance inversion method based on the sparsity constraint of the non-convex arctangent function ζ according to claim 2, wherein:
the step 2 comprises the following steps:
step 2.1: obtaining initial logarithmic wave impedance L by an initial model based on wave impedance according to a logarithmic operation formula 0 (ii) a The logarithmic operation formula is shown as formula (1):
l = ln (Z) formula (1)
Wherein Z represents wave impedance, L represents logarithmic wave impedance, and ln represents logarithmic operation;
step 2.2: constructing a non-convex arctangent function ζ as shown in equation (2):
Figure FDA0003709821330000021
wherein x represents data to be constrained, a represents the non-convex degree of an arctangent function zeta, | | represents an absolute value, and arctan represents arctangent;
step 2.3: constructing a seismic wave impedance forward model based on the sparse constraint of the non-convex arctangent function zeta, as shown in a formula (3), a formula (4) and a formula (5):
Figure FDA0003709821330000022
Figure FDA0003709821330000023
Figure FDA0003709821330000024
wherein λ represents a sparse constraint weight coefficient, β represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, | | | | | survival F Represents Frobenius norm, m represents the total sampling point of log-wave impedance, L n Represents the value of the nth sample point, w, in the log wave impedance L 1 1 st sample point, w, representing a seismic wavelet q Representing the qth sample of the seismic wavelet, q representing the seismic wavelet ensembleNumber of samples, L 0 Representing the initial log-wave impedance, a representing the degree of non-convexity of the arctangent function ζ, S 0 Representing a post-stack seismic record;
step 2.4: introducing a Lagrange multiplier term R and a dual term C on the basis of a forward model to obtain a seismic wave impedance inversion target function based on the non-convex arctangent function zeta sparse constraint, as shown in a formula (6):
Figure FDA0003709821330000025
wherein gamma is a dual term weight coefficient, gamma/lambda represents gamma divided by lambda, lambda represents a sparse constraint weight coefficient, beta represents an initial model weight coefficient, D represents a difference matrix, W represents a seismic wavelet convolution matrix, and | | Y F Represents Frobenius norm, m represents the total sampling point of log-wave impedance, L represents log-wave impedance, L 0 Representing the initial log-wave impedance, a representing the degree of non-convexity of the arctangent function ζ, S 0 Representing post-stack seismic records, R n Representing the nth value in the lagrange multiplier term R.
4. The seismic wave impedance inversion method based on the sparsity constraint of the non-convex arctangent function ζ according to claim 3, wherein:
the step 3 comprises the following steps:
step 3.1: impedance L of initial log wave 0 Initial model L as a sparsity constraint for a non-convex arctangent function ζ 1 And introducing an initial value R of a Lagrange multiplier term 1 =0 and its dual initial value C 1 =0;
Step 3.2: updating wave impedance logarithm L by using alternative direction multiplier algorithm and criterion of zero gradient extremum i+1 The calculation formula is shown in formula (7):
L i+1 =(D T W T WD+βE+γD T D) -1 [D T W T S 0 +βL 0 +γD T (R i -C i )]formula (7)
Wherein E represents an identity matrix, D T Representing the transpose of the difference matrix D, W T Representing the transpose, L, of the wavelet convolution matrix W i+1 The result of the (i + 1) th iteration, R, representing the log-wave impedance L i Representing the result of the i-th iteration of the Lagrange multiplier term R, C i Representing the ith iteration result of the even term C, T representing the transposition of the matrix, i representing the ith iteration, gamma being the weight coefficient of the even term, beta representing the weight coefficient of the initial model, L 0 Represents the initial log wave impedance, S 0 Representing a post-stack seismic record;
step 3.3: updating Lagrange multiplier term R according to alternative direction multiplier algorithm and high-order approximation method i+1 The calculation formula is shown in formula (8):
Figure FDA0003709821330000031
wherein the symbol/represents a dot division operation,
Figure FDA0003709821330000032
representing a matrix dot product operation, a representing the degree of non-convexity of an arctangent function ζ, D representing a difference matrix, C i Denotes the result of the i-th iteration, L, on the even term C i+1 The result of the (i + 1) th iteration, R, representing the log-wave impedance L i +1 Representing the (i + 1) th iteration result of the Lagrange multiplier term R;
sign represents a sign function, and the expression is shown as formula (9):
Figure FDA0003709821330000033
wherein x is j The expression vector DL i+1 +C i The jth element of (1);
f 3 and f 2 For intermediate variables, the expressions are shown in equation (10) and equation (11):
f 3 =(sqrt(f 1 .^2-f 2 .^3)-f 1 ) A (1/3) formula (10)
f 2 =b.^2./(9.*a 2 )-b./(3.*a 2 ) Formula (11)
Wherein, the symbol ^ represents a point power operation, sqrt represents a square root function, and a represents the non-convexity degree of an arctangent function ζ;
f 1 the expression of b is shown in equation (12) and equation (13):
Figure FDA0003709821330000041
b=1-a*abs(R i ) Formula (13)
Wherein R is i Representing the ith iteration result of the Lagrange multiplier term R, wherein a represents the non-convex degree of the arctangent function zeta; gamma is a dual term weight coefficient, lambda represents a sparse constraint weight coefficient, and abs represents an absolute value function;
step 3.4: updating the even term C according to the alternative direction multiplier algorithm and the Fermat lemma i+1 The calculation formula is shown as formula (14):
C i+1 =C i +(DL i+1 -R i+1 ) Formula (14)
Wherein, C i+1 Represents the i +1 th iteration result of even term C, D represents the difference matrix, C i Denotes the result of the i-th iteration, L, on the even term C i+1 Represents the (i + 1) th iteration result of the log-wave impedance L; r i+1 Represents the i +1 th iteration result of the Lagrangian multiplier term R.
5. The seismic wave impedance inversion method based on the sparsity constraint of the non-convex arctangent function ζ according to claim 4, wherein: in the step 4, whether | L is satisfied is judged i+1 -L i || 2 /||L i || 2 If yes, returning to the step 3.2 for circulation; if not, calculating the wave impedance Z according to the exponential operation r As shown in equation (15):
Figure FDA0003709821330000042
wherein the content of the first and second substances,
Figure FDA0003709821330000043
denotes the natural logarithm e as the base number, L i+1 Is an index, L i+1 Represents the i +1 th iteration result of the log-wave impedance L.
CN202210719263.2A 2022-06-23 2022-06-23 Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint Active CN114994757B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210719263.2A CN114994757B (en) 2022-06-23 2022-06-23 Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210719263.2A CN114994757B (en) 2022-06-23 2022-06-23 Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint

Publications (2)

Publication Number Publication Date
CN114994757A CN114994757A (en) 2022-09-02
CN114994757B true CN114994757B (en) 2022-12-16

Family

ID=83036990

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210719263.2A Active CN114994757B (en) 2022-06-23 2022-06-23 Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint

Country Status (1)

Country Link
CN (1) CN114994757B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115494547B (en) * 2022-10-21 2023-04-28 成都理工大学 Seismic wave impedance inversion method and system based on logarithmic total variation sparse constraint

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5798982A (en) * 1996-04-29 1998-08-25 The Trustees Of Columbia University In The City Of New York Method for inverting reflection trace data from 3-D and 4-D seismic surveys and identifying subsurface fluid and pathways in and among hydrocarbon reservoirs based on impedance models
CN110618453A (en) * 2019-08-07 2019-12-27 成都理工大学 Wave impedance inversion method based on improved damping least square method
CN112630835A (en) * 2020-12-03 2021-04-09 重庆三峡学院 High-resolution post-stack seismic wave impedance inversion method
CN113640871A (en) * 2021-08-10 2021-11-12 成都理工大学 Seismic wave impedance inversion method based on heavily-weighted L1 norm sparse constraint

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11693113B2 (en) * 2017-09-01 2023-07-04 The Trustees Of Princeton University Quantitative ultrasound imaging based on seismic full waveform inversion

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5798982A (en) * 1996-04-29 1998-08-25 The Trustees Of Columbia University In The City Of New York Method for inverting reflection trace data from 3-D and 4-D seismic surveys and identifying subsurface fluid and pathways in and among hydrocarbon reservoirs based on impedance models
CN110618453A (en) * 2019-08-07 2019-12-27 成都理工大学 Wave impedance inversion method based on improved damping least square method
CN112630835A (en) * 2020-12-03 2021-04-09 重庆三峡学院 High-resolution post-stack seismic wave impedance inversion method
CN113640871A (en) * 2021-08-10 2021-11-12 成都理工大学 Seismic wave impedance inversion method based on heavily-weighted L1 norm sparse constraint

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于稀疏约束贝叶斯估计的相对波阻抗反演;邸海滨等;《石油物探》;20110325(第02期);全文 *

Also Published As

Publication number Publication date
CN114994757A (en) 2022-09-02

Similar Documents

Publication Publication Date Title
US8700370B2 (en) Method, system and program storage device for history matching and forecasting of hydrocarbon-bearing reservoirs utilizing proxies for likelihood functions
US20210356623A1 (en) Rock Reservoir Structure Characterization Method, Device, Computer-Readable Storage Medium and Electronic Equipment
Ratto et al. State dependent parameter metamodelling and sensitivity analysis
CN112989708B (en) Well logging lithology identification method and system based on LSTM neural network
CN111368247B (en) Sparse representation regularization prestack AVO inversion method based on fast orthogonal dictionary
CN110208862B (en) Seismic inversion method based on mixed high-order fractional-order ATpV sparse regularization
CN111048163B (en) Shale oil hydrocarbon retention amount (S1) evaluation method based on high-order neural network
CN114994757B (en) Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint
CN113640871B (en) Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint
WO2013158873A2 (en) System and method for calibrating permeability for use in reservoir modeling
CN111027882A (en) Method for evaluating brittleness index by utilizing conventional logging data based on high-order neural network
CN113093272A (en) Time domain full waveform inversion method based on convolutional coding
CN112796738A (en) Stratum permeability calculation method combining array acoustic logging and conventional logging
CN110927791A (en) Method and device for predicting fluid by utilizing seismic data based on deep learning
CN114897047B (en) Multi-sensor data drift detection method based on depth dictionary
CN113009564B (en) Seismic data processing method and device
CN116090283A (en) Aviation electromagnetic three-dimensional inversion method based on compressed sensing and preconditioned random gradient
CN115494547B (en) Seismic wave impedance inversion method and system based on logarithmic total variation sparse constraint
CN108387928B (en) Data construction method based on seismic feature transformation space
CN113707216A (en) Infiltration immune cell proportion counting method
CN114755744A (en) Total organic carbon well logging interpretation method and system based on mud shale heterogeneity characteristics
CN111594156A (en) Method and system for calculating saturation of natural gas hydrate
Fu Caianiello neural network method for geophysical inverse problems
CN111435170B (en) Method and system for calculating thin layer offset profile
CN114063169B (en) Wave impedance inversion method, system, equipment and storage medium

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant