CN113640871B - Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint - Google Patents
Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint Download PDFInfo
- Publication number
- CN113640871B CN113640871B CN202110912265.9A CN202110912265A CN113640871B CN 113640871 B CN113640871 B CN 113640871B CN 202110912265 A CN202110912265 A CN 202110912265A CN 113640871 B CN113640871 B CN 113640871B
- Authority
- CN
- China
- Prior art keywords
- wave impedance
- logarithm
- updated
- seismic
- matrix
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 36
- 239000011159 matrix material Substances 0.000 claims abstract description 66
- 230000009977 dual effect Effects 0.000 claims abstract description 25
- 238000001914 filtration Methods 0.000 claims description 4
- 238000004364 calculation method Methods 0.000 description 10
- 238000010586 diagram Methods 0.000 description 5
- 238000005070 sampling Methods 0.000 description 4
- 230000000694 effects Effects 0.000 description 3
- 238000005457 optimization Methods 0.000 description 3
- 230000008602 contraction Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- OVSKIKFHRZPJSS-UHFFFAOYSA-N 2,4-D Chemical compound OC(=O)COC1=CC=C(Cl)C=C1Cl OVSKIKFHRZPJSS-UHFFFAOYSA-N 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000011478 gradient descent method Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/282—Application of seismic models, synthetic seismograms
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/30—Analysis
- G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/364—Seismic filtering
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/20—Trace signal pre-filtering to select, remove or transform specific events or signal components, i.e. trace-in/trace-out
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/63—Seismic attributes, e.g. amplitude, polarity, instant phase
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/66—Subsurface modeling
- G01V2210/665—Subsurface modeling using geostatistical modeling
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention discloses a seismic wave impedance inversion method based on a re-weighted L1 norm sparse constraint, which comprises the following steps: acquiring seismic data, and obtaining an initial wave impedance logarithm by utilizing the seismic data; respectively constructing a re-weighting matrix and a convolution matrix, and constructing a forward model based on re-weighting L1 norm sparse constraint by utilizing the initial wave impedance logarithm, the re-weighting matrix and the convolution matrix; combining Lagrange multiplier terms and dual terms on the basis of the forward model based on the re-weighted L1 norm sparse constraint to obtain a seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint; solving the seismic wave impedance inversion objective function based on the re-weighting L1 norm sparse constraint by combining an alternate direction multiplier algorithm to obtain a wave impedance logarithm, and updating the wave impedance logarithm to obtain an updated wave impedance logarithm; and obtaining an inversion result of the wave impedance by using the wave impedance logarithm and the updated wave impedance logarithm.
Description
Technical Field
The invention relates to the technical field of oil and gas exploration and seismic inversion, in particular to a seismic wave impedance inversion method based on a heavy weighting L1 norm sparse constraint.
Background
Inversion of seismic wave impedance is an important means of oil and gas exploration by obtaining the seismic wave impedance from the seismic record through an optimization method according to the mathematical and physical relationship between the post-stack seismic record and the wave impedance. Seismic wave impedance inversion technology based on sparse constraint is an important method for seismic wave impedance inversion. The sparse constraint is used as priori information to be introduced into the forward model, so that the accuracy of the 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. According to the method, an L1 norm is used as a constraint term to construct a forward model, an optimization problem formed by the forward model is subjected to iterative solution through an optimization method, and finally a wave impedance inversion result is obtained. Tayler et al used the L1 norm for seismic inversion for the first time and demonstrated that this method can improve inversion accuracy. Wang et al used L1 norm constraint for single pass wave impedance inversion and solved using gradient descent method, demonstrating the feasibility of this method in seismic wave impedance inversion. However, the L1 norm cannot fully mine sparsity, which results in a pseudo layer in the inversion result, and thus the accuracy is reduced.
Disclosure of Invention
The technical problem solved by the scheme provided by the embodiment of the invention is that the L1 norm cannot fully excavate sparsity, so that a pseudo layer appears in an inversion result, and the accuracy of the inversion result is reduced.
According to the seismic wave impedance inversion method based on the heavy weighting L1 norm sparse constraint provided by the embodiment of the invention, the method comprises the following steps:
acquiring seismic data, and obtaining an initial wave impedance logarithm by utilizing the seismic data;
respectively constructing a re-weighting matrix and a convolution matrix, and constructing a forward model based on re-weighting L1 norm sparse constraint by utilizing the initial wave impedance logarithm, the re-weighting matrix and the convolution matrix;
combining Lagrange multiplier terms and dual terms on the basis of the forward model based on the re-weighted L1 norm sparse constraint to obtain a seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint;
solving the seismic wave impedance inversion objective function based on the re-weighting L1 norm sparse constraint by combining an alternate direction multiplier algorithm to obtain a wave impedance logarithm, and updating the wave impedance logarithm to obtain an updated wave impedance logarithm;
and obtaining an inversion result of the wave impedance by using the wave impedance logarithm and the updated wave impedance logarithm.
According to the scheme provided by the embodiment of the invention, the re-weighting matrix is added on the basis of the L1 norm sparse constraint, the pseudo layer phenomenon caused by the L1 norm sparse constraint is successfully eliminated by introducing re-weighting, and the effect of improving the inversion result accuracy is achieved.
Drawings
The accompanying drawings, which are included to provide a further understanding of the invention and are incorporated in and constitute a part of this specification, illustrate embodiments of the invention and together with the description serve to explain the invention. In the drawings:
FIG. 1 is a flow chart of a seismic wave impedance inversion method based on a re-weighted L1 norm sparse constraint provided by the invention;
FIG. 2 is a block flow diagram of a seismic wave impedance inversion method based on a re-weighted L1 norm sparse constraint provided by the invention;
FIG. 3 is a noisy seismic section provided by the present invention;
FIG. 4 is a schematic diagram of an initial model of multi-channel wave impedance provided by the present invention;
FIG. 5 is a schematic diagram of the inversion result of multi-channel wave impedance provided by the invention;
FIG. 6 is a schematic diagram showing the result of updating the impedance of the single-channel side channel wave according to the present invention.
Detailed Description
The following detailed description of the preferred embodiments of the present invention is provided in conjunction with the accompanying drawings, and it is to be understood that the preferred embodiments described below are merely illustrative and explanatory of the invention, and are not restrictive of the invention.
FIG. 1 is a flowchart of a seismic wave impedance inversion method based on a weighted L1 norm sparse constraint according to an embodiment of the present invention, as shown in FIG. 1, including:
step S101: acquiring seismic data, and obtaining an initial wave impedance logarithm by utilizing the seismic data;
step S102: respectively constructing a re-weighting matrix and a convolution matrix, and constructing a forward model based on re-weighting L1 norm sparse constraint by utilizing the initial wave impedance logarithm, the re-weighting matrix and the convolution matrix;
step S103: combining Lagrange multiplier terms and dual terms on the basis of the forward model based on the re-weighted L1 norm sparse constraint to obtain a seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint;
step S104: solving the seismic wave impedance inversion objective function based on the re-weighting L1 norm sparse constraint by combining an alternate direction multiplier algorithm to obtain a wave impedance logarithm, and updating the wave impedance logarithm to obtain an updated wave impedance logarithm;
step S105: and obtaining an inversion result of the wave impedance by using the wave impedance logarithm and the updated wave impedance logarithm.
Specifically, the acquiring the seismic data and obtaining the initial wave impedance logarithm using the seismic data includes: acquiring a seismic record S containing post-stack 0 Seismic data, wavelet data w and log data; from the post-stack seismic record S 0 Extracting the position information used for representing the position of the stratum, and carrying out interpolation filtering processing on the logging data by utilizing the position information to obtain an initial model Z of wave impedance 0 The method comprises the steps of carrying out a first treatment on the surface of the Initial model Z based on wave impedance 0 Obtaining initial wave impedance logarithm L according to a logarithm operation formula 0 。
Further, the forward model based on the re-weighted L1 norm sparse constraint includes:
wherein ,the value of L is expressed as J 1 (L) minimum value, W represents a convolution matrix, D represents a difference matrix, L represents the log wave impedance to be inverted, M represents a re-weighting matrix and, I 2 The two norms of the vector are represented, I 1 Represents L1 norm, μ represents a heavily weighted variation weight coefficient, α represents an initial model Z 0 Weight coefficient of (c) in the above-mentioned formula (c).
Further, the obtaining the seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint includes:
wherein R represents a lagrangian multiplier, C represents a dual term, and λ represents a dual term weight coefficient.
Further, the logarithm of wave impedance includes:
L 1 ＝(D T W T WD+αE+λD T M 0 T M 0 D) -1 [D T W T S 0 +αL 0 +λD T M 0 T (R 1 -C 1 )]
wherein E represents an identity matrix, M 0 Represents the initial value of the re-weighting matrix, R 1 Represents the initial value of the Lagrangian multiplier, C 1 Representing the initial values of the dual terms.
Further, the updating the wave impedance logarithm to obtain an updated wave impedance logarithm includes:
using the logarithm of wave impedance L 1 Initial value R of Lagrangian multiplier term 1 Updating to obtain updated Lagrange multiplier sub-term R 2 ；
Wherein the updated Lagrangian multiplier term R 2 Comprising the following steps:
wherein ,representing a matrix point multiplication operation, sign represents a sign function, and the expression is as follows:
wherein ,xi Representing vector M 0 DL 1 +C 1 Is the i-th element of (c).
Further, the updating the wave impedance logarithm to obtain an updated wave impedance logarithm further includes:
using the logarithm of wave impedance L 1 And the updated Lagrangian multiplier term R 2 Initial value C of pair-pair item 1 Updating to obtain updated dual item C 2 ；
Wherein the updated dual item C 2 Comprising the following steps:
C 2 ＝C 1 +(M 0 DL 1 -R 2 )。
further, the updating the wave impedance logarithm to obtain an updated wave impedance logarithm further includes:
using the logarithm of wave impedance L 1 For the initial value M of the re-weighting matrix 0 Performing update processing to obtain an updated re-weighting matrix M 1 ；
Using the updated re-weighting matrix M 1 The updated Lagrangian multiplier term R 2 Said updated pair-items C 2 Logarithm of wave impedance L 1 Obtaining updated wave impedance logarithm L 2 ；
Wherein the updated re-weighting matrix M 1 Comprising the following steps:
wherein ,represents M 1 Diagonal jth element,>representation DL 1 The j-th element, ε, represents a fraction close to 0;
the updated wave impedance logarithm L 2 Comprising the following steps:
L 2 ＝(D T W T WD+αE+λD T M 1T M 1 D) -1 [D T W T S 0 +αL 0 +λD T M 1T (R 2 -C 2 )]。
further, the obtaining the inversion result of the wave impedance by using the wave impedance logarithm and the updated wave impedance logarithm includes:
determining the updated logarithm of wave impedance L 2 Whether or not to satisfy L 2 -L 1 || 2 /||L 1 || 2 ＞tol；
When the judgment result is not satisfied, the updated wave impedance logarithm L is compared with the reference value 2 Performing index operation to obtain an inversion result of the wave impedance;
where tol represents the error.
Further, the inversion result of the wave impedance includes:
FIG. 2 is a block flow diagram of the seismic wave impedance inversion method based on the re-weighted L1 norm sparse constraint, as shown in FIG. 2, comprising the following steps:
step 1: preprocessing the seismic data and obtaining an initial model of wave impedance;
step 2: calculating initial log wave impedance L from processed seismic data 0 A re-weighting matrix M is built, a forward model based on re-weighting L1 norm sparse constraint is built, and Langerhans multiplier items are introduced on the basis of the forward model to build an objective function;
step 3: constructing a re-weighting matrix initial model M 0 And is connected with L 0 As an iteration initial value, inputting a Lagrangian multiplier initial value R 1 Initial value C of dual item 1 The method comprises the steps of carrying out a first treatment on the surface of the Solving the objective function by combining the alternate direction multiplier method to obtain the wave impedance logarithm L 1 Through L 1 Updating the Lagrangian multiplier to obtain an updated Lagrangian multiplier R 2 Through L 1 ，R 2 Updating the dual item to obtain updated dual item C 2 Through L 1 Updating the re-weighting matrix to obtain an updated re-weighting matrix M 1 Through M 1 、R 2 、C 2 、L 1 Obtaining updated wave impedance logarithm L 2 。
Step 4: judging whether the value before and after updating satisfies L 2 -L 1 || 2 /||L 1 || 2 > tol, if so, then L is employed 2 and M1 As an iteration initial value, R 2 、C 2 And (3) returning to the step (3) as an initial value for circulation, and if not, carrying out exponential operation to obtain an inversion result of the wave impedance.
Preferably, the step 1 includes the steps of:
step 1.1: inputting post-stack seismic records S 0 Wavelet data w and logging data;
step 1.2: extracting horizon information for representing the position of stratum from post-stack seismic records, and performing interpolation filtering on logging data to obtain an initial model Z of wave impedance 0 。
Preferably, the step 2 includes the steps of:
step 2.1: initial model Z based on wave impedance 0 Obtaining initial logarithmic wave impedance L according to logarithmic operation formula 0 The method comprises the steps of carrying out a first treatment on the surface of the The operation formula is shown in formula 1:
L 0 ＝ln(Z 0 ) Formula (1)
Where ln represents a logarithmic operation;
step 2.2: constructing a re-weighting matrix M as shown in formula 2:
wherein ,mi The i-th element on the diagonal of the re-weighting matrix M is represented, and n represents the number of sampling points of the single-channel wave impedance.
Step 2.3: inputting a re-weighted variation weight coefficient mu, an initial model weight coefficient alpha and a difference matrix D, constructing a convolution matrix W through wavelet data W, and constructing a seismic wave impedance forward model based on re-weighted L1 sparse constraint, wherein the forward model is shown in a formula 3:
wherein ||||1 Representing the L1 norm, the differential matrix D, the convolution matrix W, L representing the log-wave impedance to be inverted,the value of L is expressed as J 1 (L) minimum; i 2 Representing a vector binary norm; the representation is as shown in equations 4, 5:
wherein ,wq The q-th sampling point of the seismic wavelet is represented, and q represents the number of wavelet sampling points.
Step 2.4: the Lagrangian multiplier term R and the dual term C are introduced on the basis of the forward model, and the seismic wave impedance inversion objective function based on the weighted L1 norm sparse constraint is obtained as shown in a formula 6:
where λ is the dual weighting coefficient.
Preferably, the step 3 includes the steps of:
step 3.1: an initial model M0 of the weighting matrix is used by the identity matrix, as shown in formula 7:
step 3.2: will initially log wave impedance L 0 As an initial model for weighting L1 norm sparse constraint, inputting an initial value M of a weighting matrix 0 And introducing the initial value R of Lagrangian multiplier term 1 =0 and its dual initial value C 1 ＝0：
Step 3.3: solving the objective function by combining the alternate direction multiplier method to obtain the wave impedance logarithm L 1 The calculation principle is shown in formula 8:
where E represents the identity matrix.
Step 3.3: updating Lagrangian multiplier term R according to alternating direction multiplier algorithm and soft threshold contraction algorithm 2 The calculation formula is shown in formula 9:
wherein ,representing a matrix point multiplication operation, sign represents a sign function, and the expression is shown in formula 10:
wherein ,xi Representing vector M 0 DL 1 +C 1 Is the i-th element of (c).
Step 3.4: updating the dual term C according to the criterion that the alternating direction multiplier algorithm and the gradient are zero extremum 2 The calculation formula is shown in formula 11:
C 2 ＝C 1 +(M 0 DL 1 -R 2 ) Formula (11)
Step 3.5: computing DL 1 Updating the re-weighting matrix M 1 The calculation expression is shown in formula (12):
wherein ,represents M 1 Diagonal jth element,>representation DL 1 The j-th element, ε, represents the fraction close to 0.
Step 3.6: according to M 1 、R 2 、C 2 、L 1 Obtaining updated wave impedance logarithm L 2 . The calculation is as shown in equation 13:
preferably, the step 4.1 judges that L is L 2 -L 1 || 2 /||L 1 || 2 > tol, if so, then L is employed 2 and M1 As an iteration initial value, R 2 、C 2 As an initial value, returning to step 3.3 for looping, ifIf not, the wave impedance is calculated according to the exponential operation, as shown in equation 14:
where tol represents the error.
Example 1
As shown in fig. 3-6, the multi-channel seismic profile wave impedance inversion is as follows:
step S1.1: inputting a plurality of post-stack seismic records S 0 Extracting horizon information through seismic records, interpolation filtering the logging data to obtain an initial model Z of wave impedance 0 。
Step S1.2: extracting single-channel seismic record s according to trace number trace 0 And single-channel wave impedance initial model z 0 。
Step S2.1: initial model z based on wave impedance 0 Obtaining initial single-channel logarithmic wave impedance l according to a logarithmic operation formula 0 The method comprises the steps of carrying out a first treatment on the surface of the The operation formula is shown in formula 31:
l 0 ＝ln(z 0 ) Formula (31)
Where ln represents a logarithmic operation;
step S2.2: a re-weighting matrix M is constructed as shown in equation 32:
wherein ,mi The i-th element on the diagonal of the re-weighting matrix M is represented, and n represents the number of sampling points of the single-channel wave impedance. Using an identity matrix as an initial value M of a weighting matrix 1 As shown in equation 33:
step S2.3: inputting a re-weighted variation weight coefficient mu, an initial model weight coefficient alpha, a differential matrix D and a convolution matrix W, and constructing a seismic wave impedance forward model based on a re-weighted L1 sparse constraint, as shown in a formula 34:
wherein ||||1 Representing the L1 norm. The differential matrix D, the convolution matrix W, is represented as shown in equations 35, 36:
wherein ,wq Representing the q-th sample point of the seismic wavelet. q represents the number of wavelet samples.
Step S2.4: the Lagrangian multiplier term R and the dual term C are introduced on the basis of the forward model, and the seismic wave impedance inversion objective function based on the weighted L1 norm sparse constraint is obtained as shown in a formula 37:
where λ is the dual weighting coefficient.
Step S3.1: will initially log the wave impedance l 0 Initial model L as a heavy weighted L1 norm sparsity constraint 1 Input the initial value M of the re-weighting matrix 1 And introducing the initial value R of Lagrangian multiplier term 1 =0 and its dual initial value C 1 ＝0；
Step S3.2: updating the wave impedance logarithm l using the criterion that the alternating direction multiplier algorithm and the gradient is zero extremum 2 The calculation principle is shown in formula 38:
l 2 ＝(D T W T WD+αE+λD T M 1T M 1 D) -1 [D T W T s 0 +αl 0 +λD T M 1T (R 1 -C 1 )]
formula (38)
Where E represents the identity matrix.
Step S3.3: updating Lagrangian multiplier term R according to alternating direction multiplier algorithm and soft threshold contraction algorithm 2 The calculation formula is shown in formula 39:
wherein ,representing a matrix dot product operation, sign represents a sign function, and the expression is shown in formula 40:
wherein ,xi Representing vector M 1 Dl 2 +C 1 Is the i-th element of (c).
Step S3.4: updating the dual term C according to the criterion that the alternating direction multiplier algorithm and the gradient are zero extremum 2 The calculation formula is shown as formula 41:
C 2 ＝C 1 +(M 1 Dl 2 -R 2 ) Formula (41)
Step S3.5: calculation of Dl 2 Updating the re-weighting matrix M 2 The calculation expression is shown in formula 42:
wherein ,represents M 2 Diagonal jth element,>representation Dl 2 The j-th element.
Step S4.1: judging I 2 -l 1 || 2 /||l 1 And if yes, returning to the step S3.1 for circulation, and if no, solving the parameter to be inverted according to the relation between the parameter to be inverted and the logarithm of the wave impedance, wherein the relation is shown in a formula 43:
step S4.2: judging whether trace is larger than the number of tracks, if not, returning to the step S1.2 for circulation, and if so, arranging all single-track records according to the number of tracks to obtain a final two-dimensional profile inversion result Z r 。
And (3) effect analysis: as shown in FIGS. 3-5, the inversion Result of the wave impedance profile can better obey the trend of the seismic record compared with the Initial model, and the correctness of the method is proved, as shown in FIG. 6, the Initial model represents the wave impedance Initial model, the Result of w_L1 represents the updated wave impedance Result, true AI represents the real wave impedance, the inversion Result is very close to the real wave impedance, no pseudo-layer phenomenon is caused, and the correctness of the method is further proved.
According to the scheme provided by the embodiment of the invention, a re-weighting matrix is introduced on the L1 norm sparse constraint, the inversion updating is carried out on the initial model by combining the alternate direction multiplier method, and the optimal inversion result is output through repeated iteration. Therefore, the false 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 is solved, and the effect of improving the inversion result accuracy is achieved.
Although the present invention has been described in detail hereinabove, the present invention is not limited thereto and various modifications may be made by those skilled in the art in accordance with the principles of the present invention. Therefore, all modifications made in accordance with the principles of the present invention should be understood as falling within the scope of the present invention.
Claims (8)
1. A seismic wave impedance inversion method based on a re-weighted L1 norm sparse constraint is characterized by comprising the following steps:
acquiring seismic data, and obtaining an initial wave impedance logarithm by utilizing the seismic data;
respectively constructing a re-weighting matrix and a convolution matrix, and constructing a forward model based on re-weighting L1 norm sparse constraint by utilizing the initial wave impedance logarithm, the re-weighting matrix and the convolution matrix;
combining Lagrange multiplier terms and dual terms on the basis of the forward model based on the re-weighted L1 norm sparse constraint to obtain a seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint;
solving the seismic wave impedance inversion objective function based on the re-weighting L1 norm sparse constraint by combining an alternate direction multiplier algorithm to obtain a wave impedance logarithm, and updating the wave impedance logarithm to obtain an updated wave impedance logarithm;
obtaining an inversion result of the wave impedance by using the wave impedance logarithm and the updated wave impedance logarithm;
wherein the forward model based on the re-weighted L1 norm sparse constraint comprises:
wherein ,the value of L is expressed as J 1 (L) minimum value, W represents a convolution matrix, D represents a difference matrix, L represents the log wave impedance to be inverted, M represents a re-weighting matrix and, I 2 The two norms of the vector are represented, I 1 Represents L1 norm, μ represents a heavily weighted variation weight coefficient, α represents an initial model Z 0 Weight coefficient of (2);
wherein the seismic wave impedance inversion objective function based on the re-weighted L1 norm sparse constraint comprises:
wherein R represents a lagrangian multiplier, C represents a dual term, and λ represents a dual term weight coefficient.
2. The method of claim 1, wherein the acquiring seismic data and using the seismic data to obtain an initial log of wave impedance comprises:
acquiring a seismic record S containing post-stack 0 Seismic data, wavelet data w and log data;
from the post-stack seismic record S 0 Extracting the position information used for representing the position of the stratum, and carrying out interpolation filtering processing on the logging data by utilizing the position information to obtain an initial model Z of wave impedance 0 ；
Initial model Z based on wave impedance 0 Obtaining initial wave impedance logarithm L according to a logarithm operation formula 0 。
3. The method of claim 2, wherein the logarithm of wave impedance comprises:
L 1 ＝(D T W T WD+αE+λD T M 0 T M 0 D) -1 [D T W T S 0 +αL 0 +λD T M 0 T (R 1 -C 1 )]
wherein E represents an identity matrix, M 0 Represents the initial value of the re-weighting matrix, R 1 Represents the initial value of the Lagrangian multiplier, C 1 Representing the initial values of the dual terms.
4. A method according to claim 3, wherein said obtaining updated logarithms of wave impedance by updating said logarithms of wave impedance comprises:
using the logarithm of wave impedance L 1 Initial value R of Lagrangian multiplier term 1 Updating to obtain updated Lagrange multiplier sub-term R 2 ；
Wherein the updated Lagrangian multiplier term R 2 Comprising the following steps:
wherein ,representing a matrix point multiplication operation, sign represents a sign function, and the expression is as follows:
wherein ,xi Representing vector M 0 DL 1 +C 1 Is the i-th element of (c).
5. The method of claim 4, wherein the obtaining the updated logarithm of the wave impedance by updating the logarithm of the wave impedance further comprises:
using the logarithm of wave impedance L 1 And the updated Lagrangian multiplier term R 2 Initial value C of pair-pair item 1 Updating to obtain updated dual item C 2 ；
Wherein the updated dual item C 2 Comprising the following steps:
C 2 ＝C 1 +(M 0 DL 1 -R 2 )。
6. the method of claim 5, wherein the obtaining the updated logarithm of the wave impedance by updating the logarithm of the wave impedance further comprises:
using the logarithm of wave impedance L 1 For the initial value M of the re-weighting matrix 0 Performing update processing to obtain updated dataRe-weighting matrix M 1 ；
Using the updated re-weighting matrix M 1 The updated Lagrangian multiplier term R 2 Said updated pair-items C 2 Logarithm of wave impedance L 1 Obtaining updated wave impedance logarithm L 2 ；
Wherein the updated re-weighting matrix M 1 Comprising the following steps:
wherein ,represents M 1 Diagonal jth element,>representation DL 1 The j-th element, ε, represents a fraction close to 0;
the updated wave impedance logarithm L 2 Comprising the following steps:
L 2 ＝(D T W T WD+αE+λD T M 1T M 1 D) -1 [D T W T S 0 +αL 0 +λD T M 1T (R 2 -C 2 )]。
7. the method of claim 6, wherein using the logarithm of wave impedance and the updated logarithm of wave impedance to obtain an inversion of wave impedance comprises:
determining the updated logarithm of wave impedance L 2 Whether or not to satisfy L 2 -L 1 || 2 /||L 1 || 2 ＞tol；
When the judgment result is not satisfied, the updated wave impedance logarithm L is compared with the reference value 2 Performing index operation to obtain an inversion result of the wave impedance;
where tol represents the error.
8. The method of claim 7, wherein the inversion of the wave impedance comprises:
。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110912265.9A CN113640871B (en) | 2021-08-10 | 2021-08-10 | Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110912265.9A CN113640871B (en) | 2021-08-10 | 2021-08-10 | Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113640871A CN113640871A (en) | 2021-11-12 |
| CN113640871B true CN113640871B (en) | 2023-09-01 |
Family
ID=78420423
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110912265.9A Active CN113640871B (en) | 2021-08-10 | 2021-08-10 | Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113640871B (en) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114994757B (en) * | 2022-06-23 | 2022-12-16 | 成都理工大学 | Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint |
| CN115494547B (en) * | 2022-10-21 | 2023-04-28 | 成都理工大学 | Seismic wave impedance inversion method and system based on logarithmic total variation sparse constraint |
Citations (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2910781A1 (en) * | 2013-04-29 | 2014-11-06 | Schlumberger Canada Limited | Deghosting with adaptive operators |
| CN108287365A (en) * | 2018-01-16 | 2018-07-17 | 中国石油大学(华东) | A kind of tri- parameter synchronization inversion methods of VSP and device based on wave equation |
| CN108535775A (en) * | 2018-03-30 | 2018-09-14 | 中国石油大学(北京) | Non-stationary seismic data sound impedance inversion method and device |
| CN109212599A (en) * | 2017-06-30 | 2019-01-15 | 中国石油化工股份有限公司 | A kind of prestack Simultaneous Retrieving method of seismic data |
| CN109298445A (en) * | 2018-09-17 | 2019-02-01 | 电子科技大学 | A kind of inverse model update method based on Gaussian Profile M-H sampling |
| CN110542923A (en) * | 2019-09-02 | 2019-12-06 | 成都理工大学 | Rapid high-precision post-stack seismic impedance inversion method |
| CN110618453A (en) * | 2019-08-07 | 2019-12-27 | 成都理工大学 | Wave impedance inversion method based on improved damping least square method |
| CN111208561A (en) * | 2020-01-07 | 2020-05-29 | 自然资源部第一海洋研究所 | Seismic acoustic wave impedance inversion method based on time-varying wavelet and curvelet transformation constraint |
| CN112147681A (en) * | 2019-06-28 | 2020-12-29 | 中国石油化工股份有限公司 | Pre-stack inversion method and system based on gamma _ Zoeppritz equation |
| CN112526599A (en) * | 2019-09-17 | 2021-03-19 | 中国石油化工股份有限公司 | Wavelet phase estimation method and system based on weighted L1 norm sparsity criterion |
| CN112578439A (en) * | 2019-09-29 | 2021-03-30 | 中国石油化工股份有限公司 | Space constraint-based seismic inversion method |
| CN112630835A (en) * | 2020-12-03 | 2021-04-09 | 重庆三峡学院 | High-resolution post-stack seismic wave impedance inversion method |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2751717C (en) * | 2008-12-22 | 2015-10-06 | Schlumberger Canada Limited | Automatic dispersion extraction of multiple time overlapped acoustic signals |
| EP3168653B1 (en) * | 2015-11-05 | 2019-07-17 | CGG Services SAS | Device and method for full waveform inversion |
-
2021
- 2021-08-10 CN CN202110912265.9A patent/CN113640871B/en active Active
Patent Citations (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2910781A1 (en) * | 2013-04-29 | 2014-11-06 | Schlumberger Canada Limited | Deghosting with adaptive operators |
| CN109212599A (en) * | 2017-06-30 | 2019-01-15 | 中国石油化工股份有限公司 | A kind of prestack Simultaneous Retrieving method of seismic data |
| CN108287365A (en) * | 2018-01-16 | 2018-07-17 | 中国石油大学(华东) | A kind of tri- parameter synchronization inversion methods of VSP and device based on wave equation |
| CN108535775A (en) * | 2018-03-30 | 2018-09-14 | 中国石油大学(北京) | Non-stationary seismic data sound impedance inversion method and device |
| CN109298445A (en) * | 2018-09-17 | 2019-02-01 | 电子科技大学 | A kind of inverse model update method based on Gaussian Profile M-H sampling |
| CN112147681A (en) * | 2019-06-28 | 2020-12-29 | 中国石油化工股份有限公司 | Pre-stack inversion method and system based on gamma _ Zoeppritz equation |
| CN110618453A (en) * | 2019-08-07 | 2019-12-27 | 成都理工大学 | Wave impedance inversion method based on improved damping least square method |
| CN110542923A (en) * | 2019-09-02 | 2019-12-06 | 成都理工大学 | Rapid high-precision post-stack seismic impedance inversion method |
| CN112526599A (en) * | 2019-09-17 | 2021-03-19 | 中国石油化工股份有限公司 | Wavelet phase estimation method and system based on weighted L1 norm sparsity criterion |
| CN112578439A (en) * | 2019-09-29 | 2021-03-30 | 中国石油化工股份有限公司 | Space constraint-based seismic inversion method |
| CN111208561A (en) * | 2020-01-07 | 2020-05-29 | 自然资源部第一海洋研究所 | Seismic acoustic wave impedance inversion method based on time-varying wavelet and curvelet transformation constraint |
| CN112630835A (en) * | 2020-12-03 | 2021-04-09 | 重庆三峡学院 | High-resolution post-stack seismic wave impedance inversion method |
Non-Patent Citations (1)
| Title |
|---|
| 杨小江.支持向量回归约束去除断层阴影研究 ——以南海珠江口盆地东部YP 油区某油田为例.《华南地震》.2021,第41卷(第41期),79-82. * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113640871A (en) | 2021-11-12 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109709603B (en) | Seismic horizon identification and tracking method and system | |
| US11693139B2 (en) | Automated seismic interpretation-guided inversion | |
| US20180247227A1 (en) | Machine learning systems and methods for data augmentation | |
| CN113640871B (en) | Seismic wave impedance inversion method based on re-weighted L1 norm sparse constraint | |
| CN112989708A (en) | Well logging lithology identification method and system based on LSTM neural network | |
| CN113687433B (en) | Bi-LSTM-based magnetotelluric signal denoising method and system | |
| CA2919543A1 (en) | System and method for estimating a reservoir parameter using joint stochastic inversion of multisource geophysical data | |
| CN114723095A (en) | Missing well logging curve prediction method and device | |
| Yang et al. | Oil logging reservoir recognition based on TCN and SA-BiLSTM deep learning method | |
| CN111722283A (en) | Stratum velocity model building method | |
| CN114966861A (en) | Seismic denoising method based on Lp pseudo-norm and gamma-norm sparse low-rank constraint | |
| Liu et al. | Lithology prediction of one-dimensional residual network based on regularization constraints | |
| CN114994757B (en) | Seismic wave impedance inversion method based on non-convex arc tangent function zeta sparse constraint | |
| Hu et al. | A hybrid CNN-LSTM machine learning model for rock mechanical parameters evaluation | |
| CN113625336A (en) | Seismic wave impedance thin layer inversion method based on full convolution neural network | |
| CN111045084B (en) | Multi-wave self-adaptive subtraction method based on prediction feature extraction | |
| Sheng et al. | Arrival-time picking of microseismic events based on MSNet | |
| CN116011338A (en) | Full waveform inversion method based on self-encoder and deep neural network | |
| Li et al. | Digital construction of geophysical well logging curves using the LSTM deep-learning network | |
| CN113642232B (en) | Intelligent inversion exploration method for surface waves, storage medium and terminal equipment | |
| CN116702577A (en) | Post-stack random seismic inversion method based on Gaussian-Markov priori constraint | |
| Jordão et al. | Using Bayesian neural networks for uncertainty assessment of ore type boundaries in complex geological models | |
| CN115494547B (en) | Seismic wave impedance inversion method and system based on logarithmic total variation sparse constraint | |
| CN110865409A (en) | Seismic wave impedance inversion method based on wave impedance low-rank regularization | |
| CN114185090B (en) | Lithology and elastic parameter synchronous inversion method and device, electronic equipment and 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 |