CN103149587B - Random coupling four-dimensional seismic inversion reservoir monitoring method and device based on grid points - Google Patents

Abstract

The invention provides a random coupling four-dimensional seismic inversion reservoir monitoring method and device based on grid points, wherein the method comprises the following steps: obtaining the variable quantity of wave impedance through two times of coupled random three-dimensional seismic inversion based on a random coupling four-dimensional seismic inversion method of a grid point; repeating the steps to obtain a plurality of wave impedance variation data volumes; obtaining a probability body of the wave impedance variation in any interval according to the plurality of wave impedance variation data bodies; and monitoring the oil reservoir change by utilizing the probability body with the wave impedance variable quantity in any interval. The invention directly uses the two-time collected post-stack seismic data as input data, thereby reducing errors caused by other operations; inverting only the variation of the wave impedance and limiting the calculated amount in a smaller range; and a probability body of the wave impedance variation in any interval is given, the variation and the uncertainty of the wave impedance are quantitatively described, and the change of the oil reservoir in the production process can be effectively reflected.

Description

Based on Random Coupling time lapse seismic inverting reservoir monitoring method and the device of net point

Technical field

The present invention relates to reservoir geophysics technology, relate in particular to the Random Coupling time lapse seismic inverting reservoir monitoring method based on net point and device.

Background technology

Stochastic seismic inversion is applicable to the area of well abundant information, can provide high-resolution inversion result and the probabilistic assessment of inverting, is a kind of important reservoir geophysics technology.Stochastic seismic inversion can be divided into the inversion method based on seismic trace and the inversion method based on net point.Time lapse seismic inversion technique utilizes twice or the seismic data of multi collect, the basis of cross equalization processing is carried out the technology that seismic inversion obtains changes in reservoir information, be that maturing field increases new reserves, develop project setting, improve the important method of recovery ratio.Time lapse seismic inverting can be divided into non-coupled inversion method, coupling inversion method and difference inversion method.Non-coupled inverting is exactly carry out separately to repeatedly seismic data the reservoir model that inverting obtains different times, and the difference then comparing them obtains the variable quantity of oil reservoir; Coupling inverting first carries out inverting to first time seismic data to obtain reservoir model, then allows this model participate in the inverting of seismic data in other in period, and finally more multiple inversion result obtains the variable quantity of oil reservoir or the variable quantity of direct acquisition oil reservoir; Difference inverting directly carries out inverting to the difference of twice seismic data, obtains the variable quantity of oil reservoir.Random time lapse seismic inverting combines the advantage of stochastic inverse and four-dimensional inverting, and many scholars have carried out research to it.

LeonBarens and AlineFranche (SEG, 2005 nd Annual Meeting collection) has carried out the random non-coupled time lapse seismic inverting numerical experiment based on road.They complete the stochastic inverse based on road to twice seismic data respectively, obtain the two cover inversion results in corresponding two periods, by contrasting the figure that crosses of this two cover inversion result wave impedance in length and breadth, point out that the method can find to produce the changes in reservoir caused when inverting uses wavelet correct.

The people such as ArildBuland propose the four-dimensional inversion method of the random difference of prestack based on seismic trace (GEOPHYSICS, 71 volume 3 phases) in 2006.The method thinks that the logarithm value of different times oil reservoir elastic parameter ratio meets Gaussian distribution, and utilize the difference of different times seismic data and the linear relationship of elastic parameter difference to build likelihood probability to distribute, and then the Posterior probability distribution that the logarithm value obtaining different times oil reservoir elastic parameter ratio meets, the expected value and standard deviation finally analytically providing above-mentioned Posterior probability distribution characterizes the uncertainty of inversion result and inversion result.

The people such as YouliQuan (SEG2008 annual meeting collection of thesis) propose the another kind of Random Coupling time lapse seismic inversion method based on seismic trace, the method, using the initial model of the result of first time inverting as second time inverting, utilizes Ensemble Kalman Filter method to obtain not oil reservoir rate pattern in the same time.

Latter two does not all consider the correlativity between inversion result Zhong Ge road based on the inverting in road.By comparison, the correlativity of each net point and surrounding net point value in three dimensions is then considered based on the method for net point.

The people such as HeleneHafslundVeire propose a kind of random time lapse seismic inversion method based on net point (GEOPHYSICS, 71 volume 5 phases) in 2006.The method belongs to difference inverting, utilize the AVO intercept of different times seismic data and the difference of gradient to build likelihood probability distribution, use Markov chain monte carlo method (MCMC) to obtain reservoir pressure, the isoparametric variable quantity of saturation degree to Posterior probability distribution sampling.Use the AVO intercept of seismic data and gradient as the input of inverting due to the method instead of directly use seismic data, the extraction error of AVO intercept and gradient may affect the precision of follow-up inverting.

T.Hong (SEG2007 annual meeting collection of thesis) etc. propose the four-dimensional inversion method of a kind of Random Coupling based on net point.The method constructs the Posterior probability distribution of the permeability comprising oil reservoir, oil saturation, pressure and a petrophysical model parameter etc., uses Markov chain Monte Carlo this Posterior probability distribution of sampling to obtain inversion result.The parameter related to due to Posterior probability distribution is more, and the method may face huge calculated amount when practical application.

LongJin (SEG2007 annual meeting collection of thesis) etc. propose the another kind of random time lapse seismic inversion method based on net point.The method, in conjunction with reservoir simulation and rock physics method, uses very fast simulated annealing algorithm to obtain the porosity model meeting time lapse seismic data and well productive frontiers.Compared with Markov chain Monte carlo algorithm, very fast simulated annealing algorithm can not the fine uncertainty that must reflect inversion result.

Realizing in process of the present invention, inventor finds that in prior art, at least there are the following problems: prior art urgently needs a kind of technical scheme compared with accurate measurements changes in reservoir.

Summary of the invention

The embodiment of the present invention provides a kind of Random Coupling time lapse seismic inverting reservoir monitoring method based on net point and device, to provide a kind of technical scheme compared with accurate measurements changes in reservoir.

On the one hand, embodiments provide a kind of Random Coupling time lapse seismic inverting reservoir monitoring method based on net point, the described Random Coupling time lapse seismic inverting reservoir monitoring method based on net point comprises:

Based on the Random Coupling time lapse seismic inversion method of net point, obtained the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling;

Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

According to described multiple wave impedance variable quantity data volume, obtain wave impedance variable quantity and be in probability volume interval arbitrarily;

Described wave impedance variable quantity is utilized to be in probability volume interval arbitrarily, monitoring changes in reservoir.

Preferably, in an embodiment of the present invention, the described Random Coupling time lapse seismic inversion method based on net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with seismic data, time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top;

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m, the Grid dimension in 3D grid Zhong Ge road is set to n,

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform;

Step 4: by the lithology value filling of Sequential Indicator Simulation by each point beyond well curve in 3D grid, put original lithology value as these;

Step 5: the wave impedance value filling being planned each point beyond well curve in 3D grid by sequential Gaussian, as the original wave impedance value of these points;

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point; If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is , new wave impedance value is ; Use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively, be denoted as syn and ; If the real seismic record that this road time 1 gathers is s ₁;

Judge whether to replace original lithology value and wave impedance value by the new lithology value of this point and wave impedance value according to Metropolis-Hasting method; Specific practice is calculated factor H ₁if, be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value; H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{syn}}}^i-s_1^i)}^2+1n\left(P_{\mathrm{litho}}\left(\mathrm{litho}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\mathrm{imp}\right)\right)]\\ -[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{syn}}^i-s_1^i)}^2+\ln \left(P_{\mathrm{litho}}\left(\stackrel{\‾}{\mathrm{litho}}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\stackrel{\‾}{\mathrm{imp}}\right)\right)]\end{array},$

Wherein, m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability, P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability;

Step 7: repeatedly repeat step 6, multiplicity is greater than 30 times, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp;

Step 8: use Reservoir simulation result, adds up the probability distribution of oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing period, is converted into the probability distribution of wave impedance variable quantity by rock physics method;

Step 9: the probability distribution of the wave impedance variable quantity using step 8 to obtain, calculates the variable quantity of its wave impedance to each net point with sequential Gaussian simulation, as the original wave impedance change value of these points;

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation; If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place; Wave impedance value and the new wave impedance variable quantity of this point obtained by step 7 again with as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂; Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity; Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1; If b is less than then still replace original value by new value, otherwise retain original value, abandon new value; H ₂be calculated as follows:

$\begin{array}{c}{H}_2=[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{\Δsyn}}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\mathrm{\Δimp}\right)\right)]\\ -[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{\Δsyn}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\stackrel{\‾}{\mathrm{\Δimp}}\right)\right)]\end{array},$

Wherein, n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeatedly repeat step 10, multiplicity is greater than 30 times, obtains a three-dimensional wave impedance variation amount data volume, thus the variable quantity passing through twice random three-dimensional seismic inversion acquisition wave impedance.

Preferably, in an embodiment of the present invention, described repetition above-mentioned steps, obtain multiple wave impedance variable quantity data volume, comprise: overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to, to the Posterior probability distribution of wave impedance variable quantity sampling N time, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.

On the other hand, embodiments provide a kind of Random Coupling time lapse seismic inverting reservoir monitoring device based on net point, the described Random Coupling time lapse seismic inverting reservoir monitoring device based on net point comprises:

Seismic inversion unit, for the Random Coupling time lapse seismic inversion method based on net point, obtains the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling; Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

Acquiring unit, for according to described multiple wave impedance variable quantity data volume, acquisition place wave impedance variable quantity is in probability volume interval arbitrarily;

Monitoring means, for utilizing described wave impedance variable quantity to be in probability volume interval arbitrarily, monitoring changes in reservoir.

Preferably, in an embodiment of the present invention, described seismic inversion unit, based on the Random Coupling time lapse seismic inversion method of net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with seismic data, time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top;

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m, the Grid dimension in 3D grid Zhong Ge road is set to n,

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform;

Step 4: by the lithology value filling of Sequential Indicator Simulation by each point beyond well curve in 3D grid, as the original lithology value of these points;

Step 5: the wave impedance value filling being planned each point beyond well curve in 3D grid by sequential Gaussian, as the original wave impedance value of these points;

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point; If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is new wave impedance value is use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively, be denoted as syn and if the real seismic record that this road time 1 gathers is s ₁;

Judge whether to replace original lithology value and wave impedance value by the new lithology value of this point and wave impedance value according to Metropolis-Hasting method; Specific practice is calculated factor H ₁if, be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value; H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{syn}}}^i-s_1^i)}^2+1n\left(P_{\mathrm{litho}}\left(\mathrm{litho}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\mathrm{imp}\right)\right)]\\ -[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{syn}}^i-s_1^i)}^2+\ln \left(P_{\mathrm{litho}}\left(\stackrel{\‾}{\mathrm{litho}}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\stackrel{\‾}{\mathrm{imp}}\right)\right)]\end{array},$

Wherein, m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability, P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability;

Step 7: repeatedly repeat step 6, multiplicity is greater than 30 times, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp;

Step 8: use Reservoir simulation result, adds up the probability distribution of oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing period, is converted into the probability distribution of wave impedance variable quantity by rock physics method;

Step 9: the probability distribution of the wave impedance variable quantity using step 8 to obtain, calculates the variable quantity of its wave impedance to each net point with sequential Gaussian simulation, as the original wave impedance change value of these points;

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation; If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place; Wave impedance value and the new wave impedance variable quantity of this point obtained by step 7 again with as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂; Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity; Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1; If b is less than , then still replace original value by new value, otherwise retain original value, abandon new value; H ₂be calculated as follows:

$\begin{array}{c}{H}_2=[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{\Δsyn}}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\mathrm{\Δimp}\right)\right)]\\ -[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{\Δsyn}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\stackrel{\‾}{\mathrm{\Δimp}}\right)\right)]\end{array},$

Wherein, n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeatedly repeat step 10, multiplicity is greater than 30 times, obtains a three-dimensional wave impedance variation amount data volume, thus the variable quantity passing through twice random three-dimensional seismic inversion acquisition wave impedance.

Preferably, in an embodiment of the present invention, described seismic inversion unit, be further used for overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to sample N time to the Posterior probability distribution of wave impedance variable quantity, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.

Technique scheme has following beneficial effect: because adopt the described Random Coupling time lapse seismic inverting reservoir monitoring method based on net point to comprise: based on the Random Coupling time lapse seismic inversion method of net point, obtained the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling; Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume; According to described multiple wave impedance variable quantity data volume, obtain wave impedance variable quantity and be in probability volume interval arbitrarily; Described wave impedance variable quantity is utilized to be in probability volume interval arbitrarily, the technological means of monitoring changes in reservoir, so reach following technique effect: compared with above-mentioned existing method, the present invention is a kind of inversion method based on net point, considers the correlativity of each net point and surrounding point value in reservoir model.Secondly, the present invention directly uses the poststack seismic data of twice collection as input data, decreases other and operates the error brought.Again, the variable quantity of the present invention's inverting wave impedance, is limited to calculated amount in less scope.Finally, the present invention uses Bayesian frame and Markov chain Monte Carlo method, provide the probability volume of wave impedance variable quantity arbitrary intervals, the variable quantity of quantitative description wave impedance and uncertainty thereof, effectively can must reflect the change of oil reservoir in production run.

Accompanying drawing explanation

In order to be illustrated more clearly in the embodiment of the present invention or technical scheme of the prior art, be briefly described to the accompanying drawing used required in embodiment or description of the prior art below, apparently, accompanying drawing in the following describes is only some embodiments of the present invention, for those of ordinary skill in the art, under the prerequisite not paying creative work, other accompanying drawing can also be obtained according to these accompanying drawings.

Fig. 1 is a kind of Random Coupling time lapse seismic inverting reservoir monitoring method flow diagram based on net point of the embodiment of the present invention;

Fig. 2 is real seismic record and the theogram schematic diagram of application example time 1 of the present invention;

Fig. 3 is real seismic record and the theogram schematic diagram of application example time 2 of the present invention;

Fig. 4 is the planimetric map that application example wave impedance variable quantity of the present invention is greater than the probability volume of 500g/cc*m/s;

Fig. 5 is the sectional view that application example wave impedance variable quantity of the present invention is greater than the probability volume of 500g/cc*m/s;

Fig. 6 is a kind of Random Coupling time lapse seismic inverting reservoir monitoring apparatus structure schematic diagram based on net point of the embodiment of the present invention.

Embodiment

Below in conjunction with the accompanying drawing in the embodiment of the present invention, be clearly and completely described the technical scheme in the embodiment of the present invention, obviously, described embodiment is only the present invention's part embodiment, instead of whole embodiments.Based on the embodiment in the present invention, those of ordinary skill in the art, not making the every other embodiment obtained under creative work prerequisite, belong to the scope of protection of the invention.

As shown in Figure 1, be a kind of Random Coupling time lapse seismic inverting reservoir monitoring method flow diagram based on net point of the embodiment of the present invention, the described Random Coupling time lapse seismic inverting reservoir monitoring method based on net point comprises:

101, based on the Random Coupling time lapse seismic inversion method of net point, the variable quantity of wave impedance is obtained by twice random three-dimensional seismic inversion of coupling;

102, repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

103, according to described multiple wave impedance variable quantity data volume, obtain wave impedance variable quantity and be in probability volume interval arbitrarily;

104, described wave impedance variable quantity is utilized to be in probability volume interval arbitrarily, monitoring changes in reservoir.

Preferably, the described Random Coupling time lapse seismic inversion method based on net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with geological data.Time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top.

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m.The Grid dimension in 3D grid Zhong Ge road is n, do like this and ensure that the temporal resolution in each road of 3D grid is all greater than Preset Time resolution.Mistiming at the bottom of the top of destination layer should not change too large within the scope of 3D grid, otherwise, should the time at the bottom of the top of suitable adjustment aim layer.

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform.

Step 4: by the statistics to well data, obtains distribution histogram and the variogram of destination layer lithology, then uses Sequential Indicator Simulation the lithology value of each point beyond well curve in 3D grid to be filled, as the original lithology value of these points.

Step 5: the histogram divide lithology to add up wave impedance to well data and variogram.The wave impedance value of each point beyond well curve in 3D grid is filled, as the original wave impedance value of these points by sequential Gaussian simulation.When simulating the wave impedance value of every bit, with reference to the lithology value of this point, the identical wave impedance value of lithology point, the wave impedance histogram of this lithology and variogram near this point is used to complete simulation.

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point.If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is new wave impedance value is use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively again, be denoted as syn and if the real seismic record that this road time 1 gathers is s ₁.Judge whether to replace original lithology value and wave impedance value by the new lithology value of this point and wave impedance value according to Metropolis-Hasting method.Specific practice is calculated factor H ₁if, be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value.H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{syn}}}^i-s_1^i)}^2+1n\left(P_{\mathrm{litho}}\left(\mathrm{litho}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\mathrm{imp}\right)\right)]\\ -[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{syn}}^i-s_1^i)}^2+\ln \left(P_{\mathrm{litho}}\left(\stackrel{\‾}{\mathrm{litho}}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\stackrel{\‾}{\mathrm{imp}}\right)\right)]\end{array}---\left(1\right)$

Wherein m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability.P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability.

Step 7: repeat 45 steps 6, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp.Step 4 to 7 achieves the three-dimensional random seismic inversion method based on net point, and this is by the Markov Chain Monte Carlo algorithm process to the Posterior probability distribution sampling of wave impedance and lithology.As shown in Figure 2, for real seismic record and the theogram schematic diagram of application example time 1 of the present invention, wherein, center section is destination layer, aterrimus is seismologic record, somber grey is just drilling the theogram of acquisition with Acoustic Impedance Data imp, and the two is very little in destination layer difference, demonstrates the reliability of this 3-d inversion.

Step 8: use Reservoir simulation result, the probability distribution of statistics oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing moment, be converted into the probability distribution of wave impedance variable quantity by rock physics method, and be a Gaussian distribution by this Probability Distribution Fitting.

Step 9: to step 1 and 2 3D grids built, after the probability distribution of the wave impedance variable quantity of known steps 8 acquisition, uses sequential Gaussian simulation each net point to be calculated to the variable quantity of its wave impedance, as the original wave impedance change value of these points.

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation.If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place.Wave impedance value and the new wave impedance variable quantity of this point obtained by step 7 again with as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂.Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity.Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1.If b is less than then still replace original value by new value, otherwise retain original value, abandon new value.H ₂factor computing formula is as follows:

$\begin{array}{c}{H}_2=[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{\Δsyn}}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\mathrm{\Δimp}\right)\right)]\\ -[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{\Δsyn}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\stackrel{\‾}{\mathrm{\Δimp}}\right)\right)]\end{array}---\left(2\right)$

Wherein n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeat 45 steps 10, obtains a three-dimensional wave impedance variation amount data volume Δ imp, thus passes through the variable quantity of twice random three-dimensional seismic inversion acquisition wave impedance.Above-mentioned steps 9 to 11 completes the process of sampling to the Posterior probability distribution of wave impedance variable quantity.The surge impedance model obtained by step 7 and the wave impedance variable quantity model that step 11 obtains are added and can obtain surge impedance model corresponding to time 2, i.e. imp+ Δ imp.As shown in Figure 3, for real seismic record and the theogram schematic diagram of application example time 2 of the present invention, wherein, center section is destination layer, aterrimus is seismologic record, somber grey is just drilling the theogram of acquisition with surge impedance model imp+ Δ imp, and the two is very little in destination layer difference, demonstrates the reliability of this four-dimensional inverting further.

Preferably, described repetition above-mentioned steps, obtain multiple wave impedance variable quantity data volume, comprise: overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to sample N time to the Posterior probability distribution of wave impedance variable quantity, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.Such as, repeat 32 above-mentioned steps 4, step 5, step 6, step 7, step 9, step 10, steps 11, obtain 32 wave impedance variable quantity bodies, then utilize these 32 wave impedance variable quantity decorum meter each point wave impedance variable quantities to be greater than any interval that 500g/cc*m/s(can get wave impedance value herein, to be greater than 500g/cc*m/s in this embodiment) probability.The planimetric map of gained probability volume and sectional view are respectively as shown in Figure 4 and Figure 5; Fig. 4 is the planimetric map that application example wave impedance variable quantity of the present invention is greater than the probability volume of 500g/cc*m/s; Fig. 5 is the sectional view (center section is destination layer) that application example wave impedance variable quantity of the present invention is greater than the probability volume of 500g/cc*m/s.

Corresponding to said method embodiment, as shown in Figure 6, be a kind of Random Coupling time lapse seismic inverting reservoir monitoring apparatus structure schematic diagram based on net point of the embodiment of the present invention, the described Random Coupling time lapse seismic inverting reservoir monitoring device based on net point comprises:

Seismic inversion unit 61, for the Random Coupling time lapse seismic inversion method based on net point, obtains the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling; Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

Acquiring unit 62, for according to described multiple wave impedance variable quantity data volume, obtains wave impedance variable quantity and is in probability volume interval arbitrarily;

Monitoring means 63, for utilizing described wave impedance variable quantity to be in probability volume interval arbitrarily, monitoring changes in reservoir.

Preferably, described seismic inversion unit 61, based on the Random Coupling time lapse seismic inversion method of net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with seismic data, time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top;

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m, the Grid dimension in 3D grid Zhong Ge road is set to n,

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform;

Step 4: by the lithology value filling of Sequential Indicator Simulation by each point beyond well curve in 3D grid, as the original lithology value of these points;

Step 5: the wave impedance value filling being planned each point beyond well curve in 3D grid by sequential Gaussian, as the original wave impedance value of these points;

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point; If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is new wave impedance value is use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively, be denoted as syn and if the real seismic record that this road time 1 gathers is s ₁;

Judge whether to replace original lithology value and wave impedance value by the new lithology value of this point and wave impedance value according to Metropolis-Hasting method; Specific practice is calculated factor H ₁if: be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value; H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{syn}}}^i-s_1^i)}^2+1n\left(P_{\mathrm{litho}}\left(\mathrm{litho}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\mathrm{imp}\right)\right)]\\ -[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{syn}}^i-s_1^i)}^2+\ln \left(P_{\mathrm{litho}}\left(\stackrel{\‾}{\mathrm{litho}}\right)\right)+\ln \left(P_{\mathrm{imp}}\left(\stackrel{\‾}{\mathrm{imp}}\right)\right)]\end{array},$

Wherein, m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability, P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability;

Step 7: repeatedly repeat step 6, multiplicity is greater than 30 times, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp;

Step 8: use Reservoir simulation result, adds up the probability distribution of oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing period, is converted into the probability distribution of wave impedance variable quantity by rock physics method;

Step 9: the probability distribution of the wave impedance variable quantity using step 8 to obtain, calculates the variable quantity of its wave impedance to each net point with sequential Gaussian simulation, as the original wave impedance change value of these points;

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation; If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place; The wave impedance value of this point obtained by step 7 again and new wave impedance variable quantity and imp+ as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂; Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity; Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1; If b is less than then still replace original value by new value, otherwise retain original value, abandon new value; H ₂be calculated as follows:

$\begin{array}{c}{H}_2=[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\stackrel{\‾}{\mathrm{\Δsyn}}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\mathrm{\Δimp}\right)\right)]\\ -[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{({\mathrm{\Δsyn}}^i-s_2^i)}^2+\ln \left(P_{\mathrm{\Δimp}}\left(\stackrel{\‾}{\mathrm{\Δimp}}\right)\right)]\end{array},$

Wherein, n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeatedly repeat step 10, multiplicity is greater than 30 times, obtains a three-dimensional wave impedance variation amount data volume, thus the variable quantity passing through twice random three-dimensional seismic inversion acquisition wave impedance.

Preferably, described seismic inversion unit 61, be further used for overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to sample N time to the Posterior probability distribution of wave impedance variable quantity, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.

Embodiment of the present invention technique scheme has following beneficial effect: because adopt the described Random Coupling time lapse seismic inverting reservoir monitoring method based on net point to comprise: based on the Random Coupling time lapse seismic inversion method of net point, obtained the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling; Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume; According to described multiple wave impedance variable quantity data volume, obtain wave impedance variable quantity and be in probability volume interval arbitrarily; Described wave impedance variable quantity is utilized to be in probability volume interval arbitrarily, the technological means of monitoring changes in reservoir, so reach following technique effect: compared with above-mentioned existing method, the present invention is a kind of inversion method based on net point, considers the correlativity of each net point and surrounding point value in reservoir model.Secondly, the present invention directly uses the poststack seismic data of twice collection as input data, decreases other and operates the error brought.Again, the variable quantity of the present invention's inverting wave impedance, is limited to calculated amount in less scope.Finally, the present invention uses Bayesian frame and Markov chain Monte Carlo method, provide the probability volume of wave impedance variable quantity arbitrary intervals, the variable quantity of quantitative description wave impedance and uncertainty thereof, effectively can must reflect the change of oil reservoir in production run.

Those skilled in the art can also recognize the various illustrative components, blocks (illustrativelogicalblock) that the embodiment of the present invention is listed, unit, and step can pass through electronic hardware, computer software, or both combinations realize.For the replaceability (interchangeability) of clear displaying hardware and software, above-mentioned various illustrative components (illustrativecomponents), unit and step have universally described their function.Such function is the designing requirement realizing depending on specific application and whole system by hardware or software.Those skilled in the art for often kind of specifically application, can use the function described in the realization of various method, but this realization can should not be understood to the scope exceeding embodiment of the present invention protection.

Various illustrative logical block described in the embodiment of the present invention, or unit can pass through general processor, digital signal processor, special IC (ASIC), field programmable gate array or other programmable logic device, discrete gate or transistor logic, discrete hardware components, or the design of above-mentioned any combination realizes or operates described function.General processor can be microprocessor, and alternatively, this general processor also can be any traditional processor, controller, microcontroller or state machine.Processor also can be realized by the combination of calculation element, such as digital signal processor and microprocessor, multi-microprocessor, and a Digital Signal Processor Core combined by one or more microprocessor, or other similar configuration any realizes.

The software module that method described in the embodiment of the present invention or the step of algorithm directly can embed hardware, processor performs or the combination of both.Software module can be stored in the storage medium of other arbitrary form in RAM storer, flash memory, ROM storer, eprom memory, eeprom memory, register, hard disk, moveable magnetic disc, CD-ROM or this area.Exemplarily, storage medium can be connected with processor, with make processor can from storage medium reading information, and write information can be deposited to storage medium.Alternatively, storage medium can also be integrated in processor.Processor and storage medium can be arranged in ASIC, and ASIC can be arranged in user terminal.Alternatively, processor and storage medium also can be arranged in the different parts in user terminal.

In one or more exemplary design, the above-mentioned functions described by the embodiment of the present invention can realize in the combination in any of hardware, software, firmware or this three.If realized in software, these functions can store on the medium with computer-readable, or are transmitted on the medium of computer-readable with one or more instruction or code form.Computer readable medium comprises computer storage medium and is convenient to make to allow computer program transfer to the telecommunication media in other place from a place.Storage medium can be that any general or special computer can the useable medium of access.Such as, such computer readable media can include but not limited to RAM, ROM, EEPROM, CD-ROM or other optical disc storage, disk storage or other magnetic storage device, or other anyly may be used for carrying or store the medium that can be read the program code of form with instruction or data structure and other by general or special computer or general or special processor.In addition, any connection can be properly termed computer readable medium, such as, if software is by a concentric cable, fiber optic cables, twisted-pair feeder, Digital Subscriber Line (DSL) or being also comprised in defined computer readable medium with wireless way for transmittings such as such as infrared, wireless and microwaves from a web-site, server or other remote resource.Described video disc (disk) and disk (disc) comprise Zip disk, radium-shine dish, CD, DVD, floppy disk and Blu-ray Disc, and disk is usually with magnetic duplication data, and video disc carries out optical reproduction data with laser usually.Above-mentioned combination also can be included in computer readable medium.

Above-described embodiment; object of the present invention, technical scheme and beneficial effect are further described; be understood that; the foregoing is only the specific embodiment of the present invention; the protection domain be not intended to limit the present invention; within the spirit and principles in the present invention all, any amendment made, equivalent replacement, improvement etc., all should be included within protection scope of the present invention.

Claims (4)

1. based on a Random Coupling time lapse seismic inverting reservoir monitoring method for net point, it is characterized in that, the described Random Coupling time lapse seismic inverting reservoir monitoring method based on net point comprises:

Based on the Random Coupling time lapse seismic inversion method of net point, obtained the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling;

Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

According to described multiple wave impedance variable quantity data volume, obtain wave impedance variable quantity and be in probability volume interval arbitrarily;

Described wave impedance variable quantity is utilized to be in probability volume interval arbitrarily, monitoring changes in reservoir;

Wherein, the described Random Coupling time lapse seismic inversion method based on net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with seismic data, time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top;

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m, the Grid dimension in 3D grid Zhong Ge road is set to n,

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform;

Step 4: by the lithology value filling of Sequential Indicator Simulation by each point beyond well curve in 3D grid, as the original lithology value of these points;

Step 5: the wave impedance value filling being planned each point beyond well curve in 3D grid by sequential Gaussian, as the original wave impedance value of these points;

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point; If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is new wave impedance value is use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively, be denoted as syn and if the real seismic record that this road time 1 gathers is s ₁;

Judge whether to replace original lithology value and wave impedance value by the new lithology value of this point and wave impedance value according to Metropolis-Hasting method; Specific practice is calculated factor H ₁if, be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value; H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=\[-\frac{1}{2n_1^n}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\stackrel{\‾}{syn}}^i-s_1^i\right)}^2+\ln \left(P_{litho}\left(litho\right)\right)+\ln \left(P_{imp}\left(imp\right)\right)\]\\ -\[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\mathrm{syn}}^i-s_1^i\right)}^2+\ln \left(P_{litho}\left(\stackrel{\‾}{litho}\right)\right)+\ln \left(P_{imp}\left(\stackrel{\‾}{imp}\right)\right)\]\end{array},$

Wherein, m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability, P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability;

Step 7: repeatedly repeat step 6, multiplicity is greater than 30 times, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp;

Step 8: use Reservoir simulation result, adds up the probability distribution of oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing period, is converted into the probability distribution of wave impedance variable quantity by rock physics method;

Step 9: the probability distribution of the wave impedance variable quantity using step 8 to obtain, calculates the variable quantity of its wave impedance to each net point with sequential Gaussian simulation, as the original wave impedance change value of these points;

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation; If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place; Wave impedance value and the new wave impedance variable quantity of this point obtained by step 7 again with as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂; Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity; Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1; If b is less than then still replace original value by new value, otherwise retain original value, abandon new value; H ₂be calculated as follows:

$\begin{array}{c}{H}_2=\[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\stackrel{\‾}{\mathrm{\Δ}syn}}^i-s_2^i\right)}^2+\ln \left(P_{\mathrm{\Δ}imp}\left(\mathrm{\Δ}imp\right)\right)\]\\ -\[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\mathrm{\Δsyn}}^i-s_2^i\right)}^2+\ln \left(P_{\mathrm{\Δ}imp}\left(\stackrel{\‾}{\mathrm{\Δ}imp}\right)\right)\]\end{array},$

Wherein, n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeatedly repeat step 10, multiplicity is greater than 30 times, obtains a three-dimensional wave impedance variation amount data volume, thus the variable quantity passing through twice random three-dimensional seismic inversion acquisition wave impedance.

2. as claimed in claim 1 based on the Random Coupling time lapse seismic inverting reservoir monitoring method of net point, it is characterized in that, repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume, comprising:

Overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to sample N time to the Posterior probability distribution of wave impedance variable quantity, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.

3. based on a Random Coupling time lapse seismic inverting reservoir monitoring device for net point, it is characterized in that, the described Random Coupling time lapse seismic inverting reservoir monitoring device based on net point comprises:

Seismic inversion unit, for the Random Coupling time lapse seismic inversion method based on net point, obtains the variable quantity of wave impedance by twice random three-dimensional seismic inversion of coupling; Repeat above-mentioned steps, obtain multiple wave impedance variable quantity data volume;

Acquiring unit, for according to described multiple wave impedance variable quantity data volume, obtains wave impedance variable quantity and is in probability volume interval arbitrarily;

Monitoring means, for utilizing described wave impedance variable quantity to be in probability volume interval arbitrarily, monitoring changes in reservoir;

Wherein, described seismic inversion unit is the Random Coupling time lapse seismic inversion method based on net point, is obtained the variable quantity of wave impedance, comprising by twice random three-dimensional seismic inversion of coupling:

Step 1: build time domain 3D grid, space sampling frequency is identical with seismic data, time of each net point of the 3D grid the superiors, time of each net point of 3D grid orlop was according to the time assignment of destination layer bottom boundary according to the time assignment at interface, destination layer top;

Step 2: Preset Time resolution ax t, pushes up the maximal value t of the mistiming end in search destination layer _m, the Grid dimension in 3D grid Zhong Ge road is set to n,

Step 3: their value to time domain, and is assigned to nearest net point by the ingrated seismic trace obtain well logging and lithology Curve transform;

Step 4: by the lithology value filling of Sequential Indicator Simulation by each point beyond well curve in 3D grid, as the original lithology value of these points;

Step 5: the wave impedance value filling being planned each point beyond well curve in 3D grid by sequential Gaussian, as the original wave impedance value of these points;

Step 6: for each net point beyond well curve, utilizes indicator Kriging to set up the lithology probability distribution P of this point _litho, utilize Sequential Indicator Simulation to obtain the new lithology value of this point, recycle the wave impedance probability distribution P that simple Ke Lijin builds this point _imp, utilize sequential Gaussian simulation to obtain the new wave impedance value of this point; If the original lithology value of this point is litho, original wave impedance value is imp, and new lithology value is new wave impedance value is use original wave impedance value and new wave impedance value to calculate the theogram in this road, place at this point respectively, be denoted as syn and if the real seismic record that this road time 1 gathers is s ₁;

Judge whether to replace original lithology value and wave impedance value by this point new lithology value and wave impedance value according to Metropolis-Hasting method; Specific practice is calculated factor H ₁if, be more than or equal to 1, then replace original value by new lithology and wave impedance value; If be less than 1, then generate a random number b, if b is less than with meeting the equally distributed tandom number generator of probability between 0 to 1 then still replace original value by new value, otherwise retain original value, abandon new value; H ₁factor computing formula is as follows:

$\begin{array}{c}{H}_1=\[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\stackrel{\‾}{syn}}^i-s_1^i\right)}^2+\ln \left(P_{litho}\left(litho\right)\right)+\ln \left(P_{imp}\left(imp\right)\right)\]\\ -\[-\frac{1}{2n_1^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\mathrm{syn}}^i-s_1^i\right)}^2+\ln \left(P_{litho}\left(\stackrel{\‾}{litho}\right)\right)+\ln \left(P_{imp}\left(\stackrel{\‾}{imp}\right)\right)\]\end{array},$

Wherein, m represents the sampling number of theogram, n ₁the noise level of the seismologic record that the expression time 1 gathers, P _litho(litho) and represent original lithology value litho and new lithology value at probability distribution P _lithoin probability, P _imp(imp) and represent original wave impedance value imp and new wave impedance value at probability distribution P _impin probability;

Step 7: repeatedly repeat step 6, multiplicity is greater than 30 times, final acquisition three-dimensional lithology data litho and Acoustic Impedance Data imp;

Step 8: use Reservoir simulation result, adds up the probability distribution of oil reservoir oil saturation and factor of porosity variable quantity between twice earthquake-capturing period, is converted into the probability distribution of wave impedance variable quantity by rock physics method;

Step 9: the probability distribution of the wave impedance variable quantity using step 8 to obtain, calculates the variable quantity of its wave impedance to each net point with sequential Gaussian simulation, as the original wave impedance change value of these points;

Step 10: to the every bit of 3D grid, utilizes simple Ke Lijin to calculate the probability distribution P of each point wave impedance variable quantity _{Δ imp}, obtain the new wave impedance variable quantity of this point with sequential Gaussian simulation; If original wave impedance variable quantity is Δ imp, new value is the wave impedance value of this point obtained by step 7 and original wave impedance variable quantity with the wave impedance value of imp+ Δ imp as this point, be the theogram Δ syn in this road, place; Wave impedance value and the new wave impedance variable quantity of this point obtained by step 7 again with as the wave impedance value of this point, then do the theogram in this road, place if the seismologic record that the time 2 gathers is s ₂; Reuse Metropolis-Hasting method to judge whether to retain new wave impedance variable quantity; Specific practice is calculated factor H ₂if, be more than or equal to 1, then replace original wave impedance variable quantity with new wave impedance variable quantity; If be less than 1, then generate random number b with the tandom number generator meeting non-uniform probability distribution between 0 to 1; If b is less than then still replace original value by new value, otherwise retain original value, abandon new value; H ₂be calculated as follows:

$\begin{array}{c}{H}_2=\[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\stackrel{\‾}{\mathrm{\Δ}syn}}^i-s_2^i\right)}^2+\ln \left(P_{\mathrm{\Δ}imp}\left(\mathrm{\Δ}imp\right)\right)\]\\ -\[-\frac{1}{2n_2^2}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{\left({\mathrm{\Δsyn}}^i-s_2^i\right)}^2+\ln \left(P_{\mathrm{\Δ}imp}\left(\stackrel{\‾}{\mathrm{\Δ}imp}\right)\right)\]\end{array},$

Wherein, n ₂the noise level of the geological data that the expression time 2 gathers, P _{Δ imp}(Δ imp) and represent that original wave impedance variable quantity and new wave impedance variable quantity are at probability distribution P respectively _{Δ imp}in probability;

Step 11: repeatedly repeat step 10, multiplicity is greater than 30 times, obtains a three-dimensional wave impedance variation amount data volume, thus the variable quantity passing through twice random three-dimensional seismic inversion acquisition wave impedance.

4. as claimed in claim 3 based on the Random Coupling time lapse seismic inverting reservoir monitoring device of net point, it is characterized in that: described seismic inversion unit, be further used for overall repetition N described step 4, step 5, step 6, step 7, step 9, step 10, step 11, be equivalent to sample N time to the Posterior probability distribution of wave impedance variable quantity, obtain N number of wave impedance variable quantity data volume, wherein N be more than or equal to 30 natural number.