CN113534262B - Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method - Google Patents
Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method Download PDFInfo
- Publication number
- CN113534262B CN113534262B CN202110704118.2A CN202110704118A CN113534262B CN 113534262 B CN113534262 B CN 113534262B CN 202110704118 A CN202110704118 A CN 202110704118A CN 113534262 B CN113534262 B CN 113534262B
- Authority
- CN
- China
- Prior art keywords
- sand
- model
- reservoir
- thickness
- earthquake
- 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 37
- 238000011161 development Methods 0.000 title claims abstract description 25
- 238000007405 data analysis Methods 0.000 title claims abstract description 17
- 238000009826 distribution Methods 0.000 claims abstract description 33
- 239000011159 matrix material Substances 0.000 claims abstract description 26
- 239000004576 sand Substances 0.000 claims abstract description 22
- 230000007704 transition Effects 0.000 claims abstract description 12
- 239000002356 single layer Substances 0.000 claims abstract description 10
- 238000004088 simulation Methods 0.000 claims abstract description 8
- 238000012512 characterization method Methods 0.000 claims abstract description 7
- 238000005516 engineering process Methods 0.000 claims abstract description 7
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 5
- 238000010606 normalization Methods 0.000 claims abstract description 4
- 238000012545 processing Methods 0.000 claims abstract description 4
- 238000012216 screening Methods 0.000 claims abstract description 4
- 238000013499 data model Methods 0.000 claims abstract description 3
- 239000010410 layer Substances 0.000 claims description 22
- 208000035126 Facies Diseases 0.000 claims description 20
- 230000008021 deposition Effects 0.000 claims description 14
- 230000001186 cumulative effect Effects 0.000 claims description 12
- 230000001902 propagating effect Effects 0.000 claims description 12
- 230000004044 response Effects 0.000 claims description 9
- 238000005070 sampling Methods 0.000 claims description 9
- 238000005315 distribution function Methods 0.000 claims description 6
- 230000014509 gene expression Effects 0.000 claims description 6
- 230000010363 phase shift Effects 0.000 claims description 6
- 239000011229 interlayer Substances 0.000 claims description 4
- 230000005540 biological transmission Effects 0.000 claims description 3
- 238000005137 deposition process Methods 0.000 claims description 3
- 238000005259 measurement Methods 0.000 claims description 3
- 238000012937 correction Methods 0.000 claims 1
- 230000015572 biosynthetic process Effects 0.000 description 10
- 238000010586 diagram Methods 0.000 description 8
- CDBYLPFSWZWCQE-UHFFFAOYSA-L Sodium Carbonate Chemical compound [Na+].[Na+].[O-]C([O-])=O CDBYLPFSWZWCQE-UHFFFAOYSA-L 0.000 description 6
- 238000010276 construction Methods 0.000 description 5
- 238000000342 Monte Carlo simulation Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 230000008901 benefit Effects 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 229910052731 fluorine Inorganic materials 0.000 description 2
- 125000001153 fluoro group Chemical group F* 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 241000446313 Lamella Species 0.000 description 1
- 244000137852 Petrea volubilis Species 0.000 description 1
- 208000027697 autoimmune lymphoproliferative syndrome due to CTLA4 haploinsuffiency Diseases 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000003595 spectral effect Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012706 support-vector machine Methods 0.000 description 1
Images
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/40—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging
- G01V1/44—Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging using generators and receivers in the same well
- G01V1/48—Processing data
- G01V1/50—Analysing data
-
- 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/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/362—Effecting static or dynamic corrections; Stacking
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Remote Sensing (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
A sand-mud interbed reservoir development zone earthquake prediction method based on big data analysis comprises the following steps: firstly, the method comprises the following steps: constructing a mass sand-mud interbedded geological model based on the sedimentary characteristics of the drilled reservoir to obtain a transition probability matrix; II, secondly, the method comprises the following steps: constructing a 0-degree incident PP wave reflection seismic record of the model in a forward algorithm and a 90-degree phase Rake wavelet simulation model; thirdly, the method comprises the following steps: screening a model trace data set highly matched with the seismic record waveform of the actual stratum through a big data analysis technology, and carrying out normalization processing; matching records in the big data model which are most similar to the actual seismic waveform according to a formula; fourthly, the method comprises the following steps: and obtaining the quantitative prediction and characterization of the earthquake of the dominant development zone of the reservoir. The invention not only realizes the earthquake forward modeling; moreover, the probability density distribution of the sand-ground ratio and the maximum single-layer sand body thickness corresponding to the corresponding geological model is statistically solved, the development degree of the reservoir is indicated, and meanwhile, the uncertainty of the interpretation is quantified.
Description
Technical Field
The invention belongs to the technical field of geophysical exploration, and particularly relates to a sand-mud interbedded reservoir development zone earthquake prediction method based on big data analysis.
Background
At present, the river facies sand shale reservoir has the geological characteristics of variable sedimentary patterns, thin sand body thickness, quick transverse change, interbedded development and the like, so that the seismic reflection characteristics of the reservoir are more complex, and the reservoir is more difficult to accurately describe. How to accurately construct the mapping relation between the sand body thickness and the seismic wave field characteristics is still the key of the thin-interbed shale reservoir prediction at the present stage.
The sand thickness (Widess) proposes a single thin layer, smaller than the tuned thickness, which can be described by a linear relationship between thickness and amplitude, as detailed in the sand thickness (Widess, 1982). Partika (Partyka) et al propose a seismic spectral decomposition method to predict the journal of lamella thickness (i.e.: partyka, 1999). Seismic attribute slicing technology, one of the key technologies of seismic sedimentology, is also an important tool for thin reservoir prediction, and has better single sand body identification effect compared with seismic profile sand body tracking, which is detailed in journal (namely: zeng,2001 li, 2015). The (Hou) and the like convert the thin interbed seismic response into a sand-ground ratio by constructing the relation between the thin interbed reservoir parameters and the seismic response and calculating the tuning factors of the thin interbed reservoir parameters and the seismic response through a support vector machine, which is detailed in the journal (namely: hou, 2019). In addition, stochastic inversion techniques such as seismic waveform inversion are also common methods for identifying thin interbed composite structures and spatial distributions using reflected waveform characteristics, which are described in the journal (e.g., chai, 2017). However, since it is difficult to cover all seismic response characteristics of complex river thin mutual reservoirs in the above method of "seismic mapping relation constructed based on limited simplified geological model", its explained ambiguity and uncertainty problem is not described effectively, and it too relies on seismic response, lacks drawbacks of geological constraints, and is not improved significantly.
Disclosure of Invention
The invention aims to provide a sand-mud interbedded reservoir development zone earthquake prediction method based on big data analysis, so as to solve the technical problems of uncertainty explained by the method, dependence on earthquake response and lack of geological constraint.
In order to achieve the purpose, the specific technical scheme of the sand-mud interbedded reservoir development zone earthquake prediction method based on big data analysis is as follows:
a method for predicting an earthquake of a sand-mud interbed reservoir development zone based on big data analysis comprises the following operation steps:
the first step is as follows: constructing a mass sand-mud interbedded geological model based on the sedimentary characteristics of the drilled reservoir;
(1) Carrying out deposition characteristic description based on a one-dimensional Markov chain model, and carrying out statistics to obtain a transition probability matrix; markov chains can model the layered sequence to capture the preferential deposition process; the basic assumptions are: the time t characteristic only depends on the underlying time t +1 characteristic, and is independent of other times; its essential characteristics can be represented by a transition probability matrix P:
(2) Endowing the geological model with thickness characteristics;
counting the thickness distribution of each sedimentary facies in a logging section corresponding to a stratum to obtain a corresponding cumulative distribution function, randomly sampling the cumulative probability to obtain a corresponding sedimentary facies thickness value, wherein the cumulative probability sampling method is represented by the formula (2):
(3) The model counts the elastic parameter distribution of each sedimentary facies in a logging section corresponding to the stratum, and the model utilizes relevant Monte Carlo to simulate and expand data volume and keeps the correlation among elastic parameters;
the second step is that: constructing a 0-degree incident PP wave reflection seismic record of the model in a propagation matrix method positive calculation method and a 90-degree phase Rake wavelet simulation model; the propagation matrix method can simulate the thickness of a thin layer to reflectionThe influence of coefficients, and the energy loss and multiple influence between layers; and in case of P-wave incidence, the sheet reflection coefficient and transmission coefficient R = [ R ] PP ，R RS ，T PP ，T PS ] T The following were used:
r＝-(A 1 -BA 2 ) -1 i P (3)
the third step: screening a model trace data set highly matched with the seismic record waveform of the actual stratum through a big data analysis technology, and carrying out normalization processing; and matching the record in the big data model which is most similar to the actual seismic waveform according to a formula (21);
the matching adopts a least square algorithm:
f(||s(d)-s(d obs )||＜ε) (21)；
the fourth step: and obtaining the quantitative prediction and characterization of the earthquake in the dominant development zone of the reservoir.
Further, in the first step, the drilled well is a typical well with a middle number in a selected work area, lithology and granularity characteristic analysis is carried out on the selected thin interbed, and a plurality of sedimentary facies are divided according to the lithology and granularity characteristic analysis, wherein the sedimentary facies are mudstone, sandy mudstone, argillaceous sandstone and black sandstone, and the granularity gradually becomes coarse from the mudstone to the sandstone.
Further, in the first step, (1) the numerical distribution of the elements in the matrix controls the directionality of the deposition, which can accurately characterize the primary deposition pattern of the formation; by continuously adjusting the initial state, massive simulation which is in accordance with the deposition mode and does not contain thickness information can be realized; and in the formula (1), P ij Column j represents the transition probability from state i to state j for row i; (2) In the corresponding sedimentary facies endowed to the model, a sedimentary mode and thickness simulation of the geological model are realized, the thickness distribution of each lithofacies in the model is ensured to accord with the actual stratum characteristics, and the accuracy of model construction is improved; finally, mass geological models conforming to the stratum sedimentary law and sedimentary facies thickness distribution can be realized;
in formula (2), i represents the ith sample; f ω () A cumulative probability distribution function representing different dephasic thicknesses; x is a radical of a fluorine atom i Represents a random number within 0 to 1; n stands for extension pointCounting;
(3) The model simulates and expands the data volume by utilizing the relevant Monte Carlo, keeps the correlation among the elastic parameters and is endowed in the model, so that the accuracy of model construction is improved; moreover, a large number of thin mutual reservoir models which accord with the geological and geophysical characteristics of the stratum are obtained.
Furthermore, in the second step, the 90-degree phase shift wavelet can realize the symmetric response of the reservoir and the waveform; the 90-degree phase shift Rake wavelet can be obtained by performing Hilbert transform on a zero-phase wavelet; and in formula (3), A 1 、A 2 And B represents the propagation matrix of the upper, lower and intermediate thin-layer media, i P Representing the incident vector.
Further, the propagation matrices of the upper and lower layer media are respectively:
in the formulas (4) and (5), omega is the circular frequency, h is the thickness of the middle thin layer, lower corner marks P and S of beta, gamma, W, Z and S correspond to a quasi-P wave and a quasi-S wave, and 1 and 2 correspond to upper and lower layers of media; the expressions of β, γ, W and Z are:
W＝p 55 (γs x +βs z ) (8)
Z＝βp 13 s x +γp 33 s z (9)
wherein, in the above formula: pv (z) 1/2 RepresentsTaking the principal square root of the complex number z, the horizontal component s of the complex slowness vector x And a vertical component s z Given by:
E＝{[(p 33 -p 55 )cos 2 θ-(p 11 -p 55 )sin 2 θ]+(p 13 +p 55 ) 2 sin 2 2θ} 1/2 (11)
wherein s is z The notation convention for the expression is:
(+, -): downward propagating pseudo-P wave, (+, +): an upwardly propagating quasi-S-wave;
(-, -): upward propagating pseudo-P wave, (, +): a downwardly propagating pseudo-S-wave;
where θ is the angle between the Z-axis and the propagation direction, p 11 、p 33 、p 55 And p 13 Is the complex modulus, both as a function of frequency;
intermediate thin layer propagation matrix:
B＝T(0)T -1 (h) (16)
in the formula
Above B is the propagation matrix for a thin layer of a single thickness h:
in the formula, h α Is a single thin inter-layer thickness and a total thin layer thickness
The incident vector is:
i P ＝iω[β P1 ，γ P1 ，-Z P1 ，-W P1 ] T (20)；
wherein in the formulas (6) and (7), rho is the density of the medium; pv represents the main value of the complex number; t (0) and T (h) respectively represent that z takes a value of 0 and h.
Further, in the third step, in formula (21), | \8230 | | | | represents a measurement distance, and epsilon represents an error parameter; and the selection of the error parameters in the formula (21) needs to be adjusted according to the matching degree, too large parameters result in too many matched invalid models, the real stratum condition cannot be represented, and too small parameters result in the matched models losing statistical significance, and the ambiguity cannot be described through statistical characteristics.
Further, the fourth step is to further obtain quantitative prediction and characterization of the earthquake of the dominant development zone of the reservoir according to the record which is obtained by matching in the third step and is most similar to the actual earthquake waveform, namely: and solving the probability density distribution of the sand-to-ground ratio and the maximum single-layer sand body thickness of the corresponding geological model, selecting the sand-to-ground ratio and the maximum single-layer sand body thickness corresponding to the maximum probability to indicate the reservoir development degree at the corresponding actual waveform position, and quantifying the uncertainty of interpretation through the probability distribution.
The sand-mud interbed reservoir development zone earthquake prediction method based on big data analysis has the following advantages:
the method combines a propagation matrix method, so that the earthquake forward modeling is realized, and the actual reservoir characterization capability of model data is improved; moreover, the probability density distribution of the sand-ground ratio and the maximum single-layer sand body thickness corresponding to the corresponding geological model is statistically solved, the development degree of the reservoir is indicated, and meanwhile, the uncertainty of the interpretation is quantified.
Drawings
FIG. 1 is a schematic flow diagram of the present invention;
FIG. 2A is a schematic illustration of a stratigraphic depositional sequence for a well 1 of the present invention; wherein, light ash represents mudstone, medium ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 2B is a schematic illustration of a stratigraphic depositional sequence for a well 2 of the present invention; wherein, light ash represents mudstone, medium ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 2C is a schematic illustration of a stratigraphic depositional sequence for a well 3 of the present invention; wherein, light ash represents mudstone, middle ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 3 is a schematic diagram of the cumulative distribution function and sampling process corresponding to the thickness of each deposition phase according to the present invention; wherein, I, II, III and IIII respectively correspond to mudstone, sand paper mudstone, argillaceous sandstone and sandstone; (which is the actual pattern on the screen)
FIG. 4A is a schematic view of an example geologic model 1 simulated by the present invention; have deposition regularity and thickness distribution characteristics similar to typical wells; wherein, light ash represents mudstone, medium ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 4B is a schematic diagram of an example geologic model 2 simulated by the present invention; have deposition regularity and thickness distribution characteristics similar to typical wells; wherein, light ash represents mudstone, medium ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 4C is a schematic diagram of an example geologic model 3 simulated by the present invention; have deposition regularity and thickness distribution characteristics similar to typical wells; wherein, light ash represents mudstone, medium ash represents sandy mudstone, deep ash represents argillaceous sandstone, and black represents sandstone; (which is the actual pattern on the screen)
FIG. 5 is a schematic diagram of the comparison between the statistics of sand-to-ground ratio of mass geological models (histogram) and the average sand-to-ground ratio of actual strata (black dashed line) according to the present invention; (which is the actual pattern on the screen)
FIG. 6A is a cross-distribution comparison of the elastic parameters (black dots) in the well and the elastic parameters (gray dots) expanded by the associated Monte Carlo simulation (taking mudstone as an example), according to the present invention, namely: a longitudinal wave velocity-transverse wave velocity diagram; (which is the actual pattern on the screen)
FIG. 6B is a cross-distribution comparison of the elastic parameters (black dots) and the expanded elastic parameters (gray dots) by the related Monte Carlo simulation (taking mudstone as an example) in the well according to the present invention, namely: longitudinal wave velocity-density. Wherein the black line represents a linear fitting relationship diagram; (which is the actual pattern on the screen)
FIG. 7A is a schematic view of a geological model according to the present invention; wherein, light-grey-mudstone, medium-grey-sandy mudstone, deep-grey-argillaceous sandstone, black-sandstone; (which is the actual pattern on the screen)
FIG. 7B is a schematic illustration of elastic distribution in terms of dephasic distribution in a geological model according to the invention; wherein, light ash-mudstone, medium ash-sandy mudstone, deep ash-argillaceous sandstone, black-sandstone; (which is the actual pattern on the screen)
FIG. 7C is a schematic illustration of a corresponding seismic record in a geological model of the invention; wherein, light ash-mudstone, medium ash-sandy mudstone, deep ash-argillaceous sandstone, black-sandstone; (which is the actual pattern on the screen)
FIG. 8A is a schematic illustration of the results of interpreting formation sand to ground ratios in accordance with the present invention; (which is the actual pattern on the screen)
FIG. 8B is a schematic illustration of an interpretation of the maximum single-layered sand thickness in a formation according to the invention; (which is the actual pattern on the screen)
FIG. 9A is a schematic representation of the formation sand to ground probability density distribution at the well 4 location of the present invention (which is the actual graph on the screen);
FIG. 9B is a schematic diagram of the probability density distribution of the maximum single sand thickness of the formation at the well 4 location of the present invention (which is an actual graph on the screen).
Detailed Description
In order to better understand the purpose, structure and function of the invention, the method for predicting the sand-mud interbed reservoir development zone earthquake based on big data analysis is described in detail below with reference to the accompanying drawings.
As shown in fig. 1-9B, the present invention employs the following steps:
the first step is as follows: constructing a mass sand-mud interbedded geological model based on the sedimentary characteristics of the drilled reservoir;
as shown in fig. 2A-2C, the well drilling is performed to select a plurality of typical wells (in this embodiment, three typical wells) in the work area, and lithology and granularity characteristics of the selected thin interbed are analyzed, and accordingly, sedimentary facies are divided, wherein light grey represents mudstone, middle grey represents sandy mudstone, dark grey represents argillaceous sandstone, black represents sandstone, and the granularity gradually becomes coarse.
(1) Accordingly, the deposition characteristic description based on the one-dimensional Markov chain model is developed, and a transition probability matrix is obtained through statistics, as shown in the following table 1:
TABLE 1 Markov chain transition probability matrix for stratigraphic correspondences
| Mudstone | Sandy mudstone | Argillaceous | Sandstone | |
| Mudstone | ||||
| 0 | 1 | 0 | 0 | |
| Sandy mudstone | 0.5 | 0 | 0.33 | 0.16 |
| |
0 | 1 | 0 | 0 |
| |
0 | 1 | 0 | 0 |
In the above table, only lithologic transition is considered and thickness characteristics are not considered, and therefore, the transition probability (diagonal element) of the same sedimentary facies is 0. It can be seen from the table that mudstone, argillaceous sandstone and sandstone can be uniformly converted into sandy mudstone, and the sandy mudstone is converted into the three types of lithofacies according to different probabilities. The numerical distribution of the elements within the matrix controls the directionality of the deposition, which can accurately characterize the primary deposition pattern of the formation. By continuously adjusting the initial state, massive simulation which is in accordance with the deposition mode and does not contain thickness information can be realized.
Markov chains may be used as a tool to simulate a layered sequence to capture a preferential deposition process. The basic assumption is that: the time t characteristic depends only on the underlying time t +1 characteristic, independent of the other times. Its essential characteristics can be represented by a transition probability matrix P:
wherein the element P ij The ith row, the jth column represents the transition probability from state i to state j.
(2) And the thickness characteristic of the geological model is given in the next step. And (3) counting the thickness distribution of each sedimentary facies in the logging section corresponding to the stratum to obtain a corresponding cumulative distribution function, and randomly sampling the cumulative probability to obtain a corresponding sedimentary facies thickness value (as shown in figure 3). And (2) endowing the corresponding sedimentary facies of the model in the step (1) with a sedimentary mode and thickness simulation of the geological model, ensuring that the thickness distribution of each lithofacies in the model accords with the actual stratum characteristics, and improving the accuracy of model construction. Finally, mass geological models which accord with the stratum sedimentary law and sedimentary facies thickness distribution can be realized;
the cumulative probability sampling method can be expressed as:
wherein i represents the ith sample; f ω () A cumulative probability distribution function representing thicknesses of different sedimentary phases; x is a radical of a fluorine atom i Represents a random number within 0 to 1; n represents the number of expansion points.
4A-4C, three simulated geological models are shown, which have similar depositional rules and thickness distribution characteristics as typical wells; carrying out sand-ground ratio statistics on the constructed mass geological model, and comparing the sand-ground ratio statistics with the sand-ground ratio in the well of the actual stratum;
as shown in fig. 5, it can be found that the geological model constructed according to the formation deposition rule and the thickness information shows better uniformity with the actual formation, and shows the advantage of geological constraint on improving the construction accuracy of the geological model.
(3) The model makes statistics of the elastic parameter distribution of each sedimentary facies in the logging section corresponding to the formation, as shown by the black dots in fig. 6A-6B. The relevant Monte Carlo is utilized to simulate and expand the data volume and keep the correlation among the elastic parameters, as shown by gray points in figures 6A-6B, and the correlation is given to the internal model (2), so that the accuracy of model construction is improved; finally, massive thin mutual reservoir models which accord with the geological and geophysical characteristics of the stratum are obtained.
The second step: forward modeling of seismic records of the model is carried out by using a propagation matrix method and 90-degree phase wavelets; constructing a 0-degree incident PP wave reflection seismic record of the model in a propagation matrix method positive calculation method and a 90-degree phase Rake wavelet simulation model; the propagation matrix method can simulate the influence of the thickness of the thin layer on the reflection coefficient, and the influence of interlayer energy loss and multiples; 90 DEG phase shift wavelets can realize symmetric response of reservoirs and waveforms. The 90-degree phase shift Rake wavelet can be obtained by performing Hilbert transform on the zero-phase wavelet.
And in case of P-wave incidence, the sheet reflection coefficient and transmission coefficient R = [ R ] PP ，R RS ，T PP ，T PS ] T The following were used:
r＝-(A 1 -BA 2 ) -1 i P (3)
wherein A is 1 、A 2 And B represents the propagation matrix of the upper, lower and intermediate thin-layer media, i P Representing the incident vector.
The propagation matrixes of the upper medium and the lower medium are respectively as follows:
in the formulas (4) and (5), ω is the circular frequency, h is the thickness of the middle thin layer, subscripts P and S of β, γ, W, Z and S correspond to pseudo-P waves and pseudo-S waves, and 1 and 2 correspond to upper and lower layers of the medium. The expressions of β, γ, W and Z are respectively:
W＝p 55 (γ sx +β sz ) (8)
Z＝βp 13 s x +γp 33 s z (9)
wherein, in the above formula: pv (z) 1/2 Representing the horizontal component s of a complex slowness vector, taking the principal value of the square root of the complex number z x And a vertical component s z Given by:
E＝{[(p 33 -p 55 )cos 2 θ-(p 11 -p 55 )sin 2 θ]+(p 13 +p 55 ) 2 sin 2 2θ} 1/2 (11)
wherein s is z The notation convention for the expression is:
(+, -): downward propagating pseudo-P wave, (+, +): an upwardly propagating quasi-S-wave;
(-, -): upward propagating pseudo-P wave, (, +): a downwardly propagating pseudo-S-wave;
where θ is the angle between the Z axis and the propagation direction, p 11 、p 33 、p 55 And p 13 Are complex moduli, all as a function of frequency.
Intermediate lamina propagation matrix:
B＝T(0)T -1 (h) (16)
in the formula
Above B is the propagation matrix for a thin layer of a single thickness h:
in the formula, h α Is a single thin interlayer thickness and a total thin layer thickness
The incident vector is:
i P ＝iω[β P1 ，γ P1 ，-Z P1 ，-W P1 ] T (20)
wherein in the formulas (6) and (7), rho is the density of the medium; pv represents the principal value of the complex number; t (0) and T (h) respectively represent that z takes a value of 0 and h.
7A-7C, FIGS. 7A-7C illustrate a geological model, elastic parameter distributions, and corresponding seismic records. The corresponding relation between the wave trough and the sand body is good, and the accuracy of the method is verified.
The third step: and screening a model trace data set highly matched with the seismic record waveform of the actual stratum by a big data analysis technology, wherein the step needs to ensure that the actual seismic record and the model seismic trace have the same sampling rate and sampling number, and normalization processing is carried out. Matching the record most similar to the actual seismic waveform according to formula (21);
the matching adopts a least square algorithm:
f(||s(d)-s(d obs )||＜ε) (21)
wherein | \8230representsthe measurement distance, and epsilon represents the error parameter.
It should be noted that the selection of the error parameter in the formula (21) needs to be adjusted according to the matching degree, too large parameter may result in too many matched invalid models, which may not represent the real formation condition, and too small parameter may result in the loss of statistical significance of the matched models, which may not describe the ambiguity through the statistical characteristics.
The fourth step: obtaining quantitative prediction and characterization of earthquake in dominant development zone of reservoir
And (4) counting probability density distribution of the sand-to-ground ratio and the maximum single-layer sand body thickness corresponding to the similar record corresponding to the geological model obtained in the third step, selecting the sand-to-ground ratio and the maximum single-layer sand body thickness corresponding to the maximum probability to indicate the reservoir development degree of the corresponding position (as shown in figures 8A-8B), indicating the reservoir development degree, and quantifying the uncertainty of explanation through the probability distribution.
FIGS. 9A-9B illustrate the probability density distributions of sand-to-land ratio and maximum single-layer sand thickness that are statistically matched to well 4 locations in FIGS. 8A-8B, enabling a description of the multi-solution problem. And the sand-ground ratio is the thickness percentage of the sandstone and the argillaceous sandstone in the stratum, and the maximum single-layer sand body thickness is the maximum thickness of the sandstone or the argillaceous sandstone in the stratum.
The above-mentioned unexplained technologies are prior art and will not be described in detail.
It is to be understood that the present invention has been described with reference to certain embodiments, and that various changes in the features and embodiments, or equivalent substitutions may be made therein by those skilled in the art without departing from the spirit and scope of the invention. In addition, many modifications may be made to adapt a particular situation or material to the teachings of the invention without departing from the essential scope thereof. Therefore, it is intended that the invention not be limited to the particular embodiment disclosed, but that the invention will include all embodiments falling within the scope of the appended claims.
Claims (5)
1. A sand-mud interbed reservoir development zone earthquake prediction method based on big data analysis is characterized by comprising the following operation steps: the first step is as follows: constructing a mass sand-mud interbedded geological model based on the sedimentary characteristics of the drilled reservoir;
(1) Carrying out deposition characteristic description based on a one-dimensional Markov chain model, and counting to obtain a transition probability matrix; the markov chain can simulate a layered sequence to capture a preferential deposition process; the basic assumption is that: the time t characteristic only depends on the underlying time t +1 characteristic, and is independent of other times; its essential characteristics can be represented by a transition probability matrix P:
(2) Endowing the geological model with thickness characteristics;
counting the thickness distribution of each sedimentary facies in a logging section corresponding to a stratum to obtain a corresponding cumulative distribution function, randomly sampling the cumulative probability to obtain a corresponding sedimentary facies thickness value, wherein the cumulative probability sampling method is represented by the formula (2):
(3) The model counts the elastic parameter distribution of each sedimentary facies in a logging section corresponding to the stratum, and the model utilizes relevant Monte Carlo to simulate and expand data volume and keeps the correlation among elastic parameters;
the second step is that: constructing a 0-degree incident PP wave reflection seismic record of the model in a propagation matrix method correction algorithm and a 90-degree phase Rake wavelet simulation model; the propagation matrix method can simulate a thin layerThe influence of thickness on reflection coefficient, and the energy loss and multiple influence between layers; and in the case of P-wave incidence, the sheet reflection coefficient and transmission coefficient R = [ R ] PP ，R RS ，T PP ，T PS ] T The following were used:
r＝-(A 1 -BA 2 ) -1 i P (3)
the third step: screening a model trace data set highly matched with the seismic record waveform of the actual stratum through a big data analysis technology, and carrying out normalization processing; and matching the record in the big data model which is most similar to the actual seismic waveform according to a formula (21);
the matching adopts a least square algorithm:
f(||s(d)-s(d obs )||＜ε) (21)；
in the third step, in a formula (21), | \8230 | | | represents a measurement distance, and epsilon represents an error parameter; the selection of the error parameters in the formula (21) needs to be adjusted according to the matching degree, too large parameters result in too many matched invalid models and the real stratum condition cannot be represented, and too small parameters result in the matched models losing statistical significance and the ambiguity cannot be described through statistical characteristics;
the fourth step: and obtaining the quantitative prediction and characterization of the earthquake in the dominant development zone of the reservoir.
2. The method for predicting the earthquake of the developmental zone of the sand-mud interbed reservoir based on the big data analysis as claimed in claim 1, wherein in the first step, the drilled wells are selected from a plurality of typical wells in a work area, lithology and granularity characteristics of the selected thin interbed are analyzed, and a plurality of sedimentary facies are divided according to the lithology and granularity characteristics, wherein the plurality of sedimentary facies are mudstones, sandy mudstones, argillaceous sandstones, and the granularity gradually becomes coarser from mudstones to sandstones.
3. The method for predicting the sand-mud interbed reservoir developmental zone earthquake according to claim 1, wherein in the second step, 90 ° phase shift rake wavelets realize symmetric response of reservoir and waveform; the 90-degree phase shift Rake wavelet can be obtained by performing Hilbert transform on the zero-phase wavelet.
4. The method for predicting the development zone earthquake of the sand-mud interbed reservoir based on the big data analysis according to claim 3, wherein propagation matrixes of upper and lower media are respectively as follows:
in the formulas (4) and (5), omega is the circular frequency, h is the thickness of the middle thin layer, lower corner marks P and S of beta, gamma, W, Z and S correspond to a quasi-P wave and a quasi-S wave, and 1 and 2 correspond to upper and lower layers of media; the expressions of β, γ, W and Z are:
W＝p 55 (γs x +βs z ) (8)
Z＝βp 13 s x +γp 33 s z (9)
wherein, in the above formula: pv (z x) 1/2 Representing the horizontal component s of the complex slowness vector, taking the principal square root of the complex number z x And a vertical component s z Given by:
E＝{[(p 33 -p 55 )cos 2 θ-(p 11 --p 55 )sin 2 θ]+(p 13 +p 55 ) 2 sin 2 2θ} 1/2 (11)
wherein s is z The notation convention for the expression is:
(+, -): downward propagating pseudo-P wave, (+, +): an upwardly propagating pseudo-S-wave;
(-, -): upward propagating pseudo-P wave, (, +): a downwardly propagating pseudo-S-wave;
where θ is the angle between the Z axis and the propagation direction, p 11 、p 33 、p 55 And p 13 Is the complex modulus, both as a function of frequency;
intermediate thin layer propagation matrix:
B＝T(0)T -1 (h) (16)
in the formula
Above B is the propagation matrix for a thin layer of a single thickness h:
in the formula, h α Is a single thin inter-layer thickness and a total thin layer thickness
The incident vector is:
i P ＝iω[β P1 ，γ P1 ，-Z P1 ，-W P1 ] T (20)。
5. the method for predicting the earthquake of the development zone of the sand-mud interbedded reservoir based on the big data analysis as claimed in claim 1, wherein the fourth step is to further obtain the quantitative prediction and characterization of the earthquake of the dominant development zone of the reservoir according to the record which is obtained by matching in the third step and is most similar to the actual earthquake waveform, namely: and solving the probability density distribution of the sand-to-ground ratio and the maximum single-layer sand body thickness of the corresponding geological model, selecting the sand-to-ground ratio and the maximum single-layer sand body thickness corresponding to the maximum probability to indicate the reservoir development degree at the corresponding actual waveform position, and quantifying the uncertainty of interpretation through the probability distribution.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110704118.2A CN113534262B (en) | 2021-06-24 | 2021-06-24 | Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110704118.2A CN113534262B (en) | 2021-06-24 | 2021-06-24 | Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113534262A CN113534262A (en) | 2021-10-22 |
| CN113534262B true CN113534262B (en) | 2023-02-03 |
Family
ID=78125789
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110704118.2A Active CN113534262B (en) | 2021-06-24 | 2021-06-24 | Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113534262B (en) |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2014169499A1 (en) * | 2013-04-19 | 2014-10-23 | 中国石油大学(华东) | Method for identifying and interpreting three-dimensional structure of ancient karst reservoir layer of carbonate rock |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| GB2397664B (en) * | 2003-01-24 | 2005-04-20 | Schlumberger Holdings | System and method for inferring geological classes |
| US8874477B2 (en) * | 2005-10-04 | 2014-10-28 | Steven Mark Hoffberg | Multifactorial optimization system and method |
| US8374836B2 (en) * | 2008-11-12 | 2013-02-12 | Geoscape Analytics, Inc. | Methods and systems for constructing and using a subterranean geomechanics model spanning local to zonal scale in complex geological environments |
| CN110618450B (en) * | 2018-06-20 | 2021-07-27 | 中国石油化工股份有限公司 | Intelligent gas-bearing property prediction method for tight reservoir based on rock physical modeling |
| CN111007567A (en) * | 2018-10-08 | 2020-04-14 | 中国石油化工股份有限公司 | Sand shale thin interbed prediction method and system based on seismic waveform inversion |
| CN111983671B (en) * | 2019-05-23 | 2023-04-07 | 中国石油天然气股份有限公司 | Shallow water lake basin reservoir prediction method and device based on micro-ancient landform restoration |
| CN112558153B (en) * | 2019-09-25 | 2022-03-29 | 中国石油天然气股份有限公司 | Oil and gas reservoir prediction method and device for two-phase medium |
| CN111487692B (en) * | 2020-04-27 | 2022-05-20 | 吉林大学 | Method for predicting seismic response characteristics and reservoir thickness of salt shale oil rhythm layer |
-
2021
- 2021-06-24 CN CN202110704118.2A patent/CN113534262B/en active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2014169499A1 (en) * | 2013-04-19 | 2014-10-23 | 中国石油大学(华东) | Method for identifying and interpreting three-dimensional structure of ancient karst reservoir layer of carbonate rock |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113534262A (en) | 2021-10-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN112505778B (en) | Three-dimensional in-situ characterization method for heterogeneity of shale storage and generation performance | |
| EP0891562B1 (en) | 3-d geologic modelling | |
| EP2111596B1 (en) | Method for generating reservoir models utilizing synthetic stratigraphic columns | |
| US8442770B2 (en) | Forming a geological model | |
| US7363163B2 (en) | Method for updating a geological reservoir model by means of dynamic data | |
| WO2019062655A1 (en) | Method and device for determining thin interlayer | |
| Dubrule | A review of stochastic models for petroleum reservoirs | |
| CN104200115B (en) | Geostatistics simulation based full-formation velocity modeling method | |
| EP3639062A1 (en) | A method for validating geological model data over corresponding original seismic data | |
| Al-Mudhafar | Statistical reservoir characterization, simulation, and optimization of field scale-gas assisted gravity drainage (GAGD) process with uncertainty assessments | |
| CN111007567A (en) | Sand shale thin interbed prediction method and system based on seismic waveform inversion | |
| CN111487692B (en) | Method for predicting seismic response characteristics and reservoir thickness of salt shale oil rhythm layer | |
| CN115877447A (en) | Reservoir prediction method for seismic restraint three-dimensional geological modeling under straight-flat combined well pattern condition | |
| CN110646850B (en) | Interlayer earthquake prediction method and device | |
| Yanshu et al. | A three-dimensional model of deep-water turbidity channel in Plutonio oilfield, Angola: From training image generation, optimization to multi-point geostatistical modelling | |
| Alves et al. | Simulation of acoustic impedance images by stochastic inversion of post-stack seismic reflection amplitudes and well data | |
| CN113534262B (en) | Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method | |
| CN109425889B (en) | Method for depicting ancient karst underground river | |
| Mehdipour et al. | The Best Scenario for Geostatistical Modeling of Porosity in the Sarvak Reservoir in an Iranian Oil Field, Using Electrofacies, Seismic Facies, and Seismic Attributes | |
| Yao et al. | spectral component geologic modeling: a new Technology for Integrating Seismic Information at the correct scale | |
| Pandey et al. | Prediction of 3-D Facies and Petrophysical Models using Seismic Inversion and Advanced Statistical Data Analytics in Midland Basin Study Area | |
| CN115032691A (en) | Probabilistic geological analysis method for sand-to-ground ratio of thin interbed in well drilling soil layer | |
| CN117345208B (en) | Quantitative characterization method and device for fracturing advantage area, electronic equipment and medium | |
| Chen | Research progress of fan delta sedimentary reservoirs for oilfield development | |
| Babasafari | New Approach to Reservoir Properties Prediction Using Petro-Elastic Inversion in a Transversely Isotropic Media |
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 |