CN113534262A - 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 PDF

Info

Publication number
CN113534262A
CN113534262A CN202110704118.2A CN202110704118A CN113534262A CN 113534262 A CN113534262 A CN 113534262A CN 202110704118 A CN202110704118 A CN 202110704118A CN 113534262 A CN113534262 A CN 113534262A
Authority
CN
China
Prior art keywords
model
sand
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.)
Granted
Application number
CN202110704118.2A
Other languages
Chinese (zh)
Other versions
CN113534262B (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.)
China National Offshore Oil Corp CNOOC
CNOOC China Ltd Tianjin Branch
Original Assignee
China National Offshore Oil Corp CNOOC
CNOOC China Ltd Tianjin Branch
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 China National Offshore Oil Corp CNOOC, CNOOC China Ltd Tianjin Branch filed Critical China National Offshore Oil Corp CNOOC
Priority to CN202110704118.2A priority Critical patent/CN113534262B/en
Publication of CN113534262A publication Critical patent/CN113534262A/en
Application granted granted Critical
Publication of CN113534262B publication Critical patent/CN113534262B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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. analysis, for interpretation, for correction
    • 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. analysis, for interpretation, for correction
    • 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

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: 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 in 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

Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method
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). The seismic attribute slicing technology, which is one of the key technologies of seismic sedimentology, is also an important tool for thin reservoir prediction, and has a better single sand body identification effect compared with seismic section 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 by using reflected waveform characteristics, which are detailed in magazines (i.e., Pei, 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 prediction 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 sand-mud interbed reservoir development zone earthquake prediction method based on big data analysis adopts 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;
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:
Figure BDA0003130504920000021
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):
Figure BDA0003130504920000022
thirdly, the elastic parameter distribution of each sedimentary facies in the logging section corresponding to the stratum is counted by the model, the data size is simulated and expanded by the model by utilizing related Monte Carlo, and the correlation among the elastic parameters is kept;
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 the influence of the thickness of the thin layer on the reflection coefficient, and the influence of interlayer energy loss and multiples; and in the case of P-wave incidence, the sheet reflection and transmission coefficients R ═ RPP,RRS,TPP,TPS]TThe following were used:
r＝-(A1-BA2)-1iP (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(dobs)||＜ε) (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 plurality of ports in a selected 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 sedimentary facies are mudstone, sandy mudstone, argillaceous sandstone and black sandstone, and the granularity gradually becomes coarse.
Furthermore, in the first step, numerical value distribution of elements in a matrix controls deposition directionality, and the method can accurately represent a main deposition mode of the stratum; 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), PijFor row i, column j represents the transition probability from state i to state j; secondly, in the corresponding sedimentary facies endowed with the model, the sedimentary mode and the 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 the model construction is improved; finally, mass geological models which accord with the stratum sedimentary rule and sedimentary facies thickness distribution can be realized;
in the formula II, i represents the ith sampling; fω() A cumulative probability distribution function representing thicknesses of different sedimentary phases; x is the number ofiRepresents a random number within 0 to 1; n represents the number of expansion points;
thirdly, the model utilizes the relevant Monte Carlo to simulate and expand the data size, maintains the correlation among the elastic parameters, and is endowed with 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 the zero-phase wavelet; and in formula (3), A1、A2And B represents the propagation matrix of the upper, lower and intermediate thin-layer media, iPRepresenting the incident vector.
Further, the propagation matrices of the upper and lower layer media are respectively:
Figure BDA0003130504920000041
Figure BDA0003130504920000042
in formulas fourth and fife, ω is the circular frequency, h is the thickness of the middle thin layer, lower corner marks P and S of β, γ, W, Z and S correspond to quasi-P waves and quasi-S waves, and 1 and 2 correspond to upper-layer and lower-layer media; the expressions of β, γ, W and Z are respectively:
Figure BDA0003130504920000043
Figure BDA0003130504920000044
W＝p55(γsx+βsz) (8)
Z＝βp13sx+γp33sz (9)
wherein, in the above formula: pv (z)1/2Representing the horizontal component s of a complex slowness vector, taking the principal value of the square root of the complex number zxAnd a vertical component szGiven by:
Figure BDA0003130504920000045
E＝{[(p33-p55)cos2θ-(p11-p55)sin2θ]+(p13+p55)2sin22θ}1/2 (11)
Figure BDA0003130504920000051
Figure BDA0003130504920000052
Figure BDA0003130504920000053
Figure BDA0003130504920000054
wherein s iszThe 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, p11、p33、p55And p13Complex moduli, all as a function of frequency;
intermediate thin layer propagation matrix:
B＝T(0)T-1(h) (16)
in the formula
Figure BDA0003130504920000055
Above B is the propagation matrix for a thin layer of a single thickness h:
Figure BDA0003130504920000056
in the formula, hαIs a single thin inter-layer thickness and a total thin layer thickness
Figure BDA0003130504920000057
The incident vector is:
iP＝iω[βP1,γP1,-ZP1,-WP1]T (20)。
further, in the third step, in formula (21), i. - | represents a measurement distance, and e 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 formation condition cannot be represented, and too small parameters result in the loss of statistical significance of the matched models, 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 the stratigraphic depositional sequence for a well 3 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. 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 fit 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 ash-mudstone, medium ash-sandy mudstone, deep ash-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 sand to land ratio in a formation according to 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 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;
as shown in fig. 2A-2C, the drilled wells are typically three wells in a selected area, lithology and granularity characteristics of the selected thin interbed are analyzed, and sedimentary facies are divided according to the lithology and granularity characteristics, wherein light ash represents mudstone, middle ash represents sandy mudstone, dark ash represents argillaceous sandstone, black represents sandstone, and the granularity gradually becomes thicker.
Carrying out the deposition characteristic description based on the one-dimensional Markov chain model, and counting to obtain a transition probability matrix, as shown in the following table 1:
TABLE 1 Markov chain transition probability matrix for stratigraphic correspondences
Mudstone Sandy mudstone Argillaceous sandstone Sandstone
Mudstone
0 1 0 0
Sandy mudstone 0.5 0 0.33 0.16
Argillaceous sandstone 0 1 0 0
Sandstone 0 1 0 0
In the above table, only the lithological transition is considered, and the thickness characteristic is 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:
Figure BDA0003130504920000091
wherein the element PijThe ith row, the jth column represents the transition probability from state i to state j.
And secondly, endowing the geological model with thickness characteristics. 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 (3) in the corresponding sedimentary facies of the model, realizing the 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 rule and sedimentary facies thickness distribution can be realized;
the cumulative probability sampling method can be expressed as:
Figure BDA0003130504920000101
wherein i represents the ith sample; fω() A cumulative probability distribution function representing thicknesses of different sedimentary phases; x is the number ofiRepresents 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.
Thirdly, the model counts the distribution of elastic parameters of each sedimentary facies in the logging section corresponding to the stratum, as shown by black points in fig. 6A-6B. The method utilizes the relevant Monte Carlo to simulate and expand the data volume and keep the correlation among the elastic parameters, as shown by gray points in FIGS. 6A-6B, and endows the gray points to the internal model, thereby not only improving the accuracy of model construction; finally, obtaining a mass thin mutual reservoir model which accords with the geological and geophysical characteristics of the stratum.
The second step is that: using a propagation matrix method forward modeling and a 90-degree phase wavelet to perform seismic record forward modeling of the model; constructing 0-degree incident PP wave reflection seismic records 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 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 the case of P-wave incidence, the sheet reflection and transmission coefficients R ═ RPP,RRS,TPP,TPS]TThe following were used:
r＝-(A1-BA2)-1iP (3)
wherein A is1、A2And B represents the propagation matrix of the upper, lower and intermediate thin-layer media, iPRepresenting the incident vector.
The propagation matrixes of the upper medium and the lower medium are respectively as follows:
Figure BDA0003130504920000111
Figure BDA0003130504920000112
in formulas fourth and fife, ω is the circular frequency, h is the thickness of the intermediate thin layer, lower corner marks P and S of β, γ, W, Z and S correspond to pseudo P-waves and pseudo S-waves, and 1 and 2 correspond to the upper and lower layer media. The expressions of β, γ, W and Z are respectively:
Figure BDA0003130504920000113
Figure BDA0003130504920000114
W＝p55(γsx+βsz) (8)
Z＝βp13sx+γp33sz (9)
wherein, in the above formula: pv (z)1/2Representing the horizontal component s of a complex slowness vector, taking the principal value of the square root of the complex number zxAnd a vertical component szGiven by:
Figure BDA0003130504920000115
E＝{[(p33-p55)cos2θ-(p11-p55)sin2θ]+(p13+p55)2sin22θ}1/2 (11)
Figure BDA0003130504920000116
Figure BDA0003130504920000121
Figure BDA0003130504920000122
Figure BDA0003130504920000123
wherein s iszThe 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, p11、p33、p55And p13Are complex moduli, all as a function of frequency.
Intermediate thin layer propagation matrix:
B＝T(0)T-1(h) (16)
in the formula
Figure BDA0003130504920000124
Above B is the propagation matrix for a thin layer of a single thickness h:
Figure BDA0003130504920000125
in the formula, hαIs a sheetThickness of thin layers, and total thickness of thin layers
Figure BDA0003130504920000126
The incident vector is:
iP＝iω[βP1,γP1,-ZP1,-WP1]T (20)
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(dobs)||＜ε) (21)
wherein, | | represents a measurement distance, and epsilon represents an 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 matched model losing statistical significance, and 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 (3) counting the probability density distribution of the sand-to-ground ratio and the maximum single-layer sand body thickness corresponding to the similarity 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 (7)

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;
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:
Figure FDA0003130504910000011
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):
Figure FDA0003130504910000012
thirdly, the elastic parameter distribution of each sedimentary facies in the logging section corresponding to the stratum is counted by the model, the data size is simulated and expanded by the model by utilizing related Monte Carlo, and the correlation among the elastic parameters is kept;
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 the influence of the thickness of the thin layer on the reflection coefficient, and the influence of interlayer energy loss and multiples; and in the case of P-wave incidence, the sheet reflection and transmission coefficients R ═ RPP,RRS,TPP,TPS]TThe following were used:
r＝-(A1-BA2)-1iP (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(dobs)||＜ε) (21)；
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 typical wells in a selected 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 sedimentary facies are mudstone, sandy mudstone, argillaceous sandstone and black sandstone, and the granularity gradually becomes coarse.
3. The method for predicting the sand-mud interbedded reservoir zone earthquake based on big data analysis according to claim 1 or 2, characterized in that in the first step, numerical distribution of elements in a matrix controls deposition directionality, which can accurately represent a main deposition pattern of the stratum; 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), PijFor row i, column j represents the transition probability from state i to state j; secondly, in the corresponding sedimentary facies endowed with the model, the sedimentary mode and the 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 the model construction is improved; finally, mass geological models which accord with the stratum sedimentary rule and sedimentary facies thickness distribution can be realized;
in the formula II, i represents the ith sampling; fω() A cumulative probability distribution function representing thicknesses of different sedimentary phases; x is the number ofiRepresents a random number within 0 to 1; n represents the number of expansion points;
thirdly, the model utilizes the relevant Monte Carlo to simulate and expand the data size, maintains the correlation among the elastic parameters, and is endowed with 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.
4. The method for predicting the sand-mud interbed reservoir developmental zone earthquake according to claim 1, wherein, in the second step, 90 ° phase shift wavelets can realize symmetric response of the reservoir and the waveform; the 90-degree phase shift Rake wavelet can be obtained by performing Hilbert transform on the zero-phase wavelet; and in formula (3), A1、A2And B represents the propagation matrix of the upper, lower and intermediate thin-layer media, iPRepresenting the incident vector.
5. The method for predicting the sand-mud interbedded reservoir developmental zone earthquake according to claim 4, wherein the propagation matrices of the upper and lower media are respectively as follows:
Figure FDA0003130504910000031
Figure FDA0003130504910000032
in formulas fourth and fife, ω is the circular frequency, h is the thickness of the middle thin layer, lower corner marks P and S of β, γ, W, Z and S correspond to quasi-P waves and quasi-S waves, and 1 and 2 correspond to upper-layer and lower-layer media; the expressions of β, γ, W and Z are respectively:
Figure FDA0003130504910000033
Figure FDA0003130504910000034
W＝p55(γsx+βsz) (8)
Z＝βp13sx+γp33sz (9)
wherein, in the above formula: pv (z)1/2Representing the horizontal component s of a complex slowness vector, taking the principal value of the square root of the complex number zxAnd a vertical component szGiven by:
Figure FDA0003130504910000035
E＝{[(p33-p55)cos2θ-(p11-p55)sin2θ]+(p13+p55)2sin22θ}1/2 (11)
Figure FDA0003130504910000041
Figure FDA0003130504910000042
Figure FDA0003130504910000043
Figure FDA0003130504910000044
wherein s iszThe 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, p11、p33、p55And p13Complex moduli, all as a function of frequency;
intermediate thin layer propagation matrix:
B＝T(0)T-1(h) (16)
in the formula
Figure FDA0003130504910000045
Above B is the propagation matrix for a thin layer of a single thickness h:
Figure FDA0003130504910000046
in the formula, hαIs a single thin inter-layer thickness and a total thin layer thickness
Figure FDA0003130504910000051
The incident vector is:
iP＝iω[βP1,γP1,-ZP1,-WP1]T (20)。
6. the method for predicting the earthquake related to the development of the sand-mud interbedded reservoir based on the big data analysis according to claim 1, wherein in the third step, in formula (21), i.e. | represents a measurement distance, and e 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 formation condition cannot be represented, and too small parameters result in the loss of statistical significance of the matched models, and the ambiguity cannot be described through statistical characteristics.
7. 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.
CN202110704118.2A 2021-06-24 2021-06-24 Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method Active CN113534262B (en)

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 true CN113534262A (en) 2021-10-22
CN113534262B 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 (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20060074825A1 (en) * 2003-01-24 2006-04-06 Piotr Mirowski System and method for inferring geological classes
US20070087756A1 (en) * 2005-10-04 2007-04-19 Hoffberg Steven M Multifactorial optimization system and method
US20100121623A1 (en) * 2008-11-12 2010-05-13 Terra Nova Sciences Llc Methods and systems for constructing and using a subterranean geomechanics model spanning local to zonal scale in complex geological environments
WO2014169499A1 (en) * 2013-04-19 2014-10-23 中国石油大学(华东) Method for identifying and interpreting three-dimensional structure of ancient karst reservoir layer of carbonate rock
CN110618450A (en) * 2018-06-20 2019-12-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
CN111487692A (en) * 2020-04-27 2020-08-04 吉林大学 Method for predicting seismic response characteristics and reservoir thickness of salt shale oil rhythm layer
CN111983671A (en) * 2019-05-23 2020-11-24 中国石油天然气股份有限公司 Shallow water lake basin reservoir prediction method and device based on micro-ancient landform restoration
CN112558153A (en) * 2019-09-25 2021-03-26 中国石油天然气股份有限公司 Oil and gas reservoir prediction method and device for two-phase medium

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20060074825A1 (en) * 2003-01-24 2006-04-06 Piotr Mirowski System and method for inferring geological classes
US20070087756A1 (en) * 2005-10-04 2007-04-19 Hoffberg Steven M Multifactorial optimization system and method
US20100121623A1 (en) * 2008-11-12 2010-05-13 Terra Nova Sciences Llc Methods and systems for constructing and using a subterranean geomechanics model spanning local to zonal scale in complex geological environments
WO2014169499A1 (en) * 2013-04-19 2014-10-23 中国石油大学(华东) Method for identifying and interpreting three-dimensional structure of ancient karst reservoir layer of carbonate rock
CN110618450A (en) * 2018-06-20 2019-12-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
CN111983671A (en) * 2019-05-23 2020-11-24 中国石油天然气股份有限公司 Shallow water lake basin reservoir prediction method and device based on micro-ancient landform restoration
CN112558153A (en) * 2019-09-25 2021-03-26 中国石油天然气股份有限公司 Oil and gas reservoir prediction method and device for two-phase medium
CN111487692A (en) * 2020-04-27 2020-08-04 吉林大学 Method for predicting seismic response characteristics and reservoir thickness of salt shale oil rhythm layer

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
LEANDRO PASSOS DE FIGUEIREDO ET AL.: "Multimodal Markov chain Monte Carlo method for nonlinear petrophysical seismic inversion", 《GEOPHYSICS》 *
LEI LIU ET AL.: "The research of seismic data processing sequence in gas cloud zone-applied in Bohai P oilfield", 《2017 CGS/SEG INTERNATIONAL GEOPHYSICAL CONFERENCE》 *
刘财等: "地震波形反演技术在砂泥岩薄互层结构表征中的应用", 《地球物理学报》 *
印兴耀等: "岩石物理驱动下地震流体识别研究", 《中国科学:地球科学》 *
夏同星等: "渤海湾 Ｘ 油田气云区地震资料关键处理技术研究", 《石油物探》 *
田玉昆等: "基于马尔科夫随机场的岩性识别方法", 《地球物理学报》 *
胡求红等: "马尔科夫链分析在东海陆架盆地花港组沉积微相分析中的应用", 《地质与资源》 *

Also Published As

Publication number Publication date
CN113534262B (en) 2023-02-03

Similar Documents

Publication Publication Date Title
AU2018340369B2 (en) Method and device for determining thin interlayer
CN112505778B (en) Three-dimensional in-situ characterization method for heterogeneity of shale storage and generation performance
EP0891562B1 (en) 3-d geologic modelling
US10067253B2 (en) Method for determining sedimentary facies using 3D seismic data
EP2111596B1 (en) Method for generating reservoir models utilizing synthetic stratigraphic columns
CN104200115B (en) Geostatistics simulation based full-formation velocity modeling method
Yasin et al. Estimation of petrophysical parameters from seismic inversion by combining particle swarm optimization and multilayer linear calculator
EP3639062A1 (en) A method for validating geological model data over corresponding original seismic data
Zahmatkesh et al. Estimating Vsand and reservoir properties from seismic attributes and acoustic impedance inversion: A case study from the Mansuri oilfield, SW Iran
Al-Mudhafar How is multiple-point geostatistics of lithofacies modeling assisting for fast history matching? A case study from a sand-rich fluvial depositional environment of Zubair formation in South Rumaila oil field
Maurya et al. Qualitative and quantitative comparison of the genetic and hybrid genetic algorithm to estimate acoustic impedance from post-stack seismic data of Blackfoot field, Canada
CN115877447A (en) Reservoir prediction method for seismic restraint three-dimensional geological modeling under straight-flat combined well pattern condition
Alves et al. Simulation of acoustic impedance images by stochastic inversion of post-stack seismic reflection amplitudes and well data
CN112764110A (en) Clustered seismic facies analysis method based on limiting Boltzmann machine feature coding
CN113534262B (en) Big data analysis-based sand-mud interbed reservoir development zone earthquake prediction method
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
CN109425889B (en) Method for depicting ancient karst underground river
Li A Bayesian approach to causal and evidential analysis for uncertainty quantification throughout the reservoir forecasting process
Yao et al. spectral component geologic modeling: a new Technology for Integrating Seismic Information at the correct scale
Sadeghi et al. Global stochastic seismic inversion using turning bands simulation and co-simulation
CN115032691A (en) Probabilistic geological analysis method for sand-to-ground ratio of thin interbed in well drilling soil layer
Chen Research progress of fan delta sedimentary reservoirs for oilfield development
Yang Holistic strategies for prediction uncertainty quantification of contaminant transport and reservoir production in field cases
Esmaeilpour Geomechanical characterization of hydrocarbon reservoirs using seismic inversion and downhole measurements
Wu et al. The Application of Stochastic Seismic Inversion Method of Geostatistical Inversion in Putaohua X Area

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