CN114879259B - Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer - Google Patents

Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer Download PDF

Info

Publication number
CN114879259B
CN114879259B CN202210398462.8A CN202210398462A CN114879259B CN 114879259 B CN114879259 B CN 114879259B CN 202210398462 A CN202210398462 A CN 202210398462A CN 114879259 B CN114879259 B CN 114879259B
Authority
CN
China
Prior art keywords
seismic
saturation
fracture
velocity
hydrate layer
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210398462.8A
Other languages
Chinese (zh)
Other versions
CN114879259A (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.)
Qingdao Marine Science And Technology Center
Institute of Oceanology of CAS
Ocean University of China
Original Assignee
Institute of Oceanology of CAS
Ocean University of China
Qingdao National Laboratory for Marine Science and Technology Development Center
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 Institute of Oceanology of CAS, Ocean University of China, Qingdao National Laboratory for Marine Science and Technology Development Center filed Critical Institute of Oceanology of CAS
Priority to CN202210398462.8A priority Critical patent/CN114879259B/en
Publication of CN114879259A publication Critical patent/CN114879259A/en
Application granted granted Critical
Publication of CN114879259B publication Critical patent/CN114879259B/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/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/62—Physical property of subsurface
    • G01V2210/622—Velocity, density or impedance
    • G01V2210/6222—Velocity; travel time
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/62—Physical property of subsurface
    • G01V2210/622—Velocity, density or impedance
    • G01V2210/6224—Density
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/62—Physical property of subsurface
    • G01V2210/622—Velocity, density or impedance
    • G01V2210/6226—Impedance
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/62—Physical property of subsurface
    • G01V2210/624—Reservoir parameters
    • G01V2210/6244—Porosity

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

The invention relates to a method for predicting anisotropic seismic saturation of a fracture filling type natural gas hydrate layer. Aiming at the characteristics of the fracture filling type hydrate layer, the seismic wave impedance data volume of the fracture filling type hydrate layer is obtained through seismic inversion of seismic data and logging data; calculating a velocity, density and porosity data volume of the hydrate layer through a density and velocity low-frequency model obtained by seismic inversion; calculating the speed of a theoretical fracture filling type hydrate layer by using a porosity data volume and an anisotropic rock physical model; and comparing the speeds of actual analysis and theoretical calculation to obtain the saturation of the fracture filling type hydrate layer. The method can solve the problem that the seismic saturation of the previous isotropic rock physical model is too large, so that the prediction precision of the hydrate layer type is improved, and technical support can be provided for drilling selection of the hydrate layer type.

Description

Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer
Technical Field
The invention relates to the technical field of marine natural gas hydrate exploration, in particular to an anisotropic seismic saturation prediction method for a fracture filling type natural gas hydrate layer.
Background
Natural gas hydrate is widely distributed in sediments of deep water basins and land frozen soil zones at the edges of continents in six occurrence forms of cementing particles, cladding particles, bearing matrix/particles, pore filling, inclusion bodies, tuberculosis/crack filling and the like; in marine sediments, however, pressure coring results show that pore filling and particle-displaced fracture filling are two common occurring morphologies of natural gas hydrates. Since the pore filling type hydrate is discovered in south China sea for the first time in 2007, the research on the hydrate in aspects of seismic survey and explanation, well logging evaluation, rock physical model, reservoir formation system and the like is very mature; in recent years, some fracture filling type stratums are newly found in south China sea, and related well logging data evaluation work is carried out, but further improvement is needed in seismic data evaluation.
When the well logging data evaluation of the fracture filling type hydrate layer is carried out, the fact that the hydrate stratum presents obvious anisotropic characteristics is found, and the hydrate saturation calculated by using a traditional isotropic evaluation method is obviously higher than a pressure coring estimation result; similarly, when seismic data is used for evaluating a fracture filling type hydrate stratum, if an isotropic method is used, the amount of hydrate resources is estimated too high, so that an anisotropic rock physical model is needed to predict and evaluate the fracture filling type hydrate stratum. Therefore, the invention provides an anisotropic seismic saturation prediction method for the fracture filling type natural gas hydrate layer.
Disclosure of Invention
Aiming at the problem that the seismic prediction saturation of the fracture filling type hydrate layer is higher than the pressure core estimation result, the invention aims to improve the seismic calculation precision of the fracture filling type hydrate layer.
According to the characteristic that the crack filling type hydrate layer presents anisotropy, the technical scheme adopted by the invention for realizing the purpose is as follows:
the anisotropic seismic saturation prediction method of the fracture filling type natural gas hydrate layer comprises the following steps:
s1, generating seismic data of synthetic records by using velocity and density logging data through a convolution method;
s2, obtaining a seismic wave impedance data volume, a velocity and density low-frequency model of the fracture filling type hydrate layer through post-stack seismic wave impedance inversion by using seismic data, velocity and density logging data, and calculating a velocity, density and porosity data volume of the hydrate layer;
s3, constructing an anisotropic rock physical model of the fracture filling type hydrate by utilizing a Backus average and particle contact model, and calculating the theoretical speed of a theoretical fracture filling type hydrate layer by utilizing a porosity data body and the anisotropic rock physical model;
and S4, comparing the inversion speed with the theoretical calculation speed to estimate the saturation of the fracture filling type hydrate.
The synthetic seismic record includes:
calculating a reflection coefficient sequence r (i) of different depths of the stratum by using stratum depth domain data of actual logging speed and density:
Figure SMS_1
wherein the parameter is depth threshold, r (i) is the reflection coefficient of the ith layer, V i+1 、V i Logging speeds, p, for the i +1 th and i-th layers, respectively i+1 、ρ i Logging densities of the i +1 th layer and the i-th layer respectively;
and convolving with the input Rake wavelet information to generate a time domain synthetic seismic record:
Figure SMS_2
in the formula, S (t) is synthetic seismic record of a time domain, d tau is a time interval, w (tau) is wavelets in offset trace concentration, and r (t-tau) is a time domain reflection coefficient sequence after velocity conversion of a depth domain reflection coefficient sequence r (i).
The post-stack seismic wave impedance inversion comprises: well logging data processing, seismic data processing, synthetic seismic recording, well seismic calibration, low-frequency model establishment, and model analysis and updating.
And calculating the speed, density and porosity data volume, wherein the calculation process comprises the following steps: according to the velocity and density low-frequency model obtained by inversion of the post-stack seismic wave impedance, dividing the velocity and density low-frequency model by the wave impedance data volume to obtain a density data volume and a velocity data volume; and calculating a porosity data volume using the density data volume:
Figure SMS_3
where phi is the porosity, rho s For particle density, p is obtained by inversionDensity, p w Is the density of water.
The step S3 includes:
calculating a Lame coefficient of a saturated water stratum by using a particle contact model by taking a porosity data volume and regional stratum lithology as input parameters;
and calculating the speed of the fracture filling type hydrate layer by using Backus average by taking the Lame coefficients of the saturated water stratum and the hydrate as input parameters.
The particle contact model is used for calculating the elastic parameters of the saturated water stratum, and the formula comprises the following steps:
bulk modulus K of dry skeleton of sediment dry And shear modulus G dry ，
Figure SMS_4
Or
Figure SMS_5
Wherein φ is porosity; phi is a c Porosity that is dense random inclusion (about 40%); k HM And G HM Needs to be obtained by Hertz-Mindlin theoretical calculation; calculating the bulk modulus K and the shear modulus G of the solid particles according to mineral components by using a Hill average formula;
bulk modulus K of saturated water deposit sat And shear modulus G sat Calculated using the Gassmann equation:
Figure SMS_6
G sat ＝G dry ；
in the formula, K w Is the bulk modulus of water;
calculating the Lame coefficient lambda and mu of the saturated water deposit according to the formula:
λ＝K sat -2μ/3,
μ＝G sat .
the Backus is averagely divided into two end members of a saturation stratum and a hydrate, and the calculation formula comprises the following steps:
Figure SMS_7
Figure SMS_8
Figure SMS_9
Figure SMS_10
N＝<μ>,
ρ＝<ρ>,
Figure SMS_11
in the formula, A, C, F, L, N, rho and Q are elastic parameters in a rigidity matrix respectively, lambda and mu are Lame coefficients of a saturated water formation and a hydrate respectively, and gamma is an included angle between a wave vector and a symmetry axis; the sharp brackets represent the weighted average expression form:
<G>≡(η 1 G 1 +η 2 G 2 )
Figure SMS_12
in the formula eta 1 And η 2 G corresponds to the calculation formula in the tip bracket in the above stiffness matrix formula for the volumes of the two end members;
calculating the longitudinal and transverse wave speeds of the crack filling type hydrate layer by using the elastic parameters,
Figure SMS_13
Figure SMS_14
Figure SMS_15
in the formula, V p 、
Figure SMS_16
And &>
Figure SMS_17
The longitudinal wave velocity, the vertical polarization transverse wave and the horizontal polarization transverse wave are respectively, namely the theoretical calculation velocity.
The estimating the saturation of the fracture-filled hydrate comprises:
and calculating the error between the inversion speed and the theoretical calculation speed under different hydrate saturation conditions, wherein when the error is less than 1e-2, the saturation is the finally calculated saturation of the fracture filling type hydrate layer.
The anisotropic seismic saturation prediction device of the fracture filling type hydrate layer predicts and estimates the anisotropic seismic saturation of the fracture filling type natural gas hydrate layer through the following program modules:
the seismic data synthesis module is used for generating synthetic recorded seismic data by utilizing the velocity and density logging data through a convolution method;
the inversion and porosity data module is used for obtaining a seismic wave impedance data volume, a velocity and density low-frequency model of the fracture filling type hydrate layer by using seismic data, velocity and density logging data and performing post-stack seismic wave impedance inversion, and calculating the velocity, density and porosity data volume of the hydrate layer;
the theoretical velocity calculation module is used for constructing an anisotropic rock physical model of the fracture filling type hydrate by utilizing the Backus average and particle contact model and calculating the theoretical velocity of a theoretical fracture filling type hydrate layer by utilizing the porosity data volume and the anisotropic rock physical model;
the saturation prediction estimation module is used for comparing the inversion speed with the theoretical calculation speed to estimate the saturation of the fracture filling type hydrate;
and the visual plotting module is used for plotting each data volume in the process of predicting and estimating the saturation and visually displaying.
The anisotropic seismic saturation prediction equipment for the fracture-filled natural gas hydrate layer comprises a processing part and a storage part, wherein the storage part stores a program, and the processing part loads the program to execute the method steps to predict and estimate the anisotropic seismic saturation of the fracture-filled natural gas hydrate layer.
Compared with the prior art, the invention has the following advantages and beneficial effects:
1. the method can solve the problem of over-high prediction of seismic saturation of the fracture filling type hydrate layer, thereby improving the prediction precision of seismic data.
2. The method can also be used for evaluating the seismic saturation of the fracture filling type hydrate layer, and can be used for calculating the anisotropic saturation of three-dimensional seismic data and predicting the spatial distribution of the hydrate.
Drawings
FIG. 1 is a general flow diagram of the present invention;
FIG. 2 is a synthetic seismic record map and simulated seismic data;
FIG. 3 is a low frequency model of wave impedance, velocity and density;
FIG. 4 is a wave impedance, density, velocity and porosity data volume for post-stack seismic inversion;
FIG. 5 is a graph of saturation prediction for a fracture-filled hydrate layer according to the invention.
Detailed Description
The present invention will be described in further detail with reference to the accompanying drawings and examples.
As shown in FIG. 1, the inventive method is divided into four steps.
Step one, a synthetic seismic record part is generated by utilizing a convolution model.
1.1, according to the stratum depth domain data of the density (figure 2 (a)) and the velocity (figure 2 (b)) obtained by logging, calculating a reflection coefficient sequence r (i) of different stratum depths by using a formula (1), and performing convolution (see a formula (2)) with given Rake wavelet information to generate a time domain synthetic seismic record, wherein the red line is shown in a figure 2 (c).
Figure SMS_18
The parameter in equation (1) is the depth threshold, r (i) is the reflection coefficient of the i-th layer, V i+1 、V i Logging speed, p, for the i +1 and i layers i+1 、ρ i Logging density for the i +1 th and ith layers, respectively.
Figure SMS_19
In the formula (2), S (t) is the synthetic seismic record of the time domain, d tau is the time interval, w (tau) is the wavelet in the offset trace set, and r (t-tau) is the time domain reflection coefficient sequence after the depth domain reflection coefficient sequence r (i) is subjected to velocity conversion.
1.2, the synthetic seismic record is used and copied to further generate the seismic data as shown in fig. 2 (e).
Step two, performing a seismic inversion part by using the velocity, the density and the synthetic seismic data S (t);
and 2.1, comparing and analyzing the synthetic seismic data S (t) (such as a black line in a figure 2 (c)) and well-side seismic trace data (such as a red line in a figure 2 (c)) at a well drilling position, and acquiring the time and depth relation of the seismic data and the logging data so as to establish a wave impedance (figure 3 (a)), a velocity (figure 3 (b)) and a density low-frequency model (figure 3 (c)) as shown in figure 3.
And 2.2, continuously updating the initial model according to the established wave impedance low-frequency model as the initial model (figure 3 (a)), performing convolution on the updated model and wavelets extracted from the seismic data to generate a synthetic seismic record, comparing the synthetic seismic record with actual input data, and obtaining a wave impedance data volume shown in figure 4 (a) by final inversion when the synthetic seismic record and the actual input data meet certain precision.
2.3, the inverted density and velocity data volumes shown in fig. 4 (b) and 4 (c) can be obtained by dividing the wave impedance data volume obtained by inversion (fig. 4 (a)) by the velocity and density low-frequency models shown in fig. 3 (b) and 3 (c), respectively. Wherein the "divide by" operation is: and dividing the x and y two-dimensional coordinates of the dividend data body and the divisor data body respectively.
2.4, calculating a porosity data volume according to the density data volume;
Figure SMS_20
where phi is the porosity, rho s Is the particle density, ρ is the density obtained by inversion, ρ w Is the density of water.
Step three, calculating the theoretical speed of the fracture filling type hydrate layer;
3.1, calculating the Lame coefficient of the saturated water stratum by using a particle contact model (formula 4-9) by taking a porosity data volume and the lithology of the regional stratum (corresponding to the volume occupied by the solid particle component) as input parameters;
bulk modulus K of dry skeleton of sediment dry And shear modulus G dry ，
Figure SMS_21
Or
Figure SMS_22
Wherein φ is porosity; phi is a c Is the porosity of the dense random inclusion (about 40%). Effective bulk modulus K of skeleton HM And shear modulus G HM It needs to be calculated by Hertz-Mindlin theory:
Figure SMS_23
Figure SMS_24
wherein ν is a mineral poisson's ratio calculated using the bulk modulus K and shear modulus G of the solid particles; n is the coordination number, here taken to be 8.5; p is the effective pressure:
P＝(1-φ)(ρ s -ρ f )gh (6)
in the formula, ρ s And ρ f Density of the solid and fluid phases, respectively; g is the acceleration of gravity; h is the depth below the seafloor.
The bulk modulus K and the shear modulus G of the solid particles are calculated from the mineral components using the Hill average formula:
Figure SMS_25
Figure SMS_26
in the formula, K i And G i Respectively the bulk modulus and shear modulus of each solid particle, f i Is the volume occupied by the corresponding solid particle fraction.
Bulk modulus K of saturated water deposit sat And shear modulus G sat Calculated using the Gassmann equation:
Figure SMS_27
G sat ＝G dry ； (8)
in the formula, K w Is the bulk modulus of water.
Calculating the Lame coefficient lambda and mu of the saturated water deposit according to the formula:
λ＝K sat -2μ/3,
μ＝G sat . (9)
3.2, calculating the speed of the fracture filling type hydrate layer by using the average backsus with the Lame coefficient of the saturated water stratum and the hydrate as an input parameter;
Figure SMS_28
Figure SMS_29
Figure SMS_30
Figure SMS_31
N＝<μ>,
rho＝<ρ>,
Figure SMS_32
in the formula, A, C, F, L, N, rho and Q are elastic parameters in a rigidity matrix respectively, lambda and mu are Lame coefficients of a saturated water stratum and a hydrate respectively, rho is the density of the saturated water stratum and the hydrate respectively, and gamma is an included angle between a wave vector and a symmetry axis. The tip brackets represent a weighted average, expressed as:
<G>≡(η 1 G 1 +η 2 G 2 )
Figure SMS_33
in the formula eta 1 And η 2 G corresponds to the calculation formula in the tip bracket in the above stiffness matrix formula for the volumes of the two end members;
calculating the longitudinal and transverse wave speeds of the crack filling type hydrate layer by using the elastic parameters, namely the theoretical speed:
Figure SMS_34
Figure SMS_35
Figure SMS_36
in the formula, V p 、
Figure SMS_37
And &>
Figure SMS_38
Longitudinal wave velocity, vertically polarized transverse wave and horizontally polarized transverse wave.
And step four, calculating the saturation of the fracture filling type hydrate layer.
4.1, calculating theoretical speeds of fracture filling type hydrate layers under different hydrate saturation conditions by continuously increasing the hydrate saturation, and calculating an error between a speed (figure 4 (c)) in an inversion speed body and the theoretical speed (a result of a formula 12), wherein when the error is less than or equal to 1e-2, the saturation at the moment is the final inversion hydrate saturation;
and 4.2, calculating the hydrate saturation of each path in the inversion velocity body (figure 4 (c)) according to the steps, and finally generating a saturation data body shown in figure 5.
As shown in fig. 2, (a) is density logging; (b) velocity logging; (c) Performing convolution on a reflection coefficient sequence obtained for the longitudinal waves and the density and the seismic wavelets (d) to obtain a comparison relation graph of synthetic seismic data (black lines) and actual seismic traces (black lines) of well positions, and obtaining the depth relation of the time of the seismic data and the logging data (namely the corresponding relation of the time in the seismic data and the depth in the logging data); (d) wavelets extracted from (e) in actual seismic traces; (e) The multi-channel post-stack seismic data generated by copying the single-channel synthetic seismic data are also seismic input data for post-stack inversion.
As shown in FIG. 3, (a), (b) and (c) are initial low frequency models of wave impedance, velocity and density, respectively, established during seismic inversion. The low frequency model is from the relationship between the time of the seismic data and the logging depth established in fig. 2 (c), namely, the logging density and velocity variation curve with the depth can be converted into a variation curve with time matched with the seismic data through the relationship, and then the low frequency model with continuous space, velocity, density and wave impedance variation with time is formed through the transverse interpolation of the curve.
As shown in fig. 4, (a) is a wave impedance data volume obtained by inversion of well logging data processing, seismic data processing, synthetic seismic recording, well seismic calibration, low frequency model building, model analysis and updating, etc.; (b) The density data volume is obtained by dividing the wave impedance data volume (fig. 4 (a)) by the velocity low frequency model (fig. 3 (b)); (c) The velocity data volume is obtained by dividing the wave impedance data volume (fig. 4 (a)) by the density low-frequency model (fig. 3 (c)); (d) Is a porosity data volume calculated from the density data volume (fig. 4 (b)) using equation (3).
As shown in fig. 5, (a) a saturation data volume for a conventional isotropic calculation; (b) Is the calculated saturation data volume for anisotropy in the present invention. The hydrate saturation estimated by actual pressure coring is about 20%, and as can be seen from the figure, the saturation highest value of the traditional isotropic calculation is 50%, which is obviously higher, and the anisotropic saturation is closer to the actual result.
The above description is only a preferred embodiment of the present invention, and is not intended to limit the present invention, and all simple modifications, changes and equivalent structural changes made to the above embodiment according to the technical spirit of the present invention still fall within the protection scope of the technical solution of the present invention.

Claims (10)

1. The anisotropic seismic saturation prediction method of the fracture filling type natural gas hydrate layer is characterized by comprising the following steps of:
s1, generating seismic data of a synthetic seismic record by using velocity and density logging data through a convolution method;
s2, obtaining a seismic wave impedance data volume, a velocity and density low-frequency model of the fracture filling type hydrate layer through post-stack seismic wave impedance inversion by using seismic data, velocity and density logging data, and calculating a velocity, density and porosity data volume of the hydrate layer;
s3, constructing an anisotropic rock physical model of the fracture filling type hydrate by utilizing a Backus average and particle contact model, and calculating the theoretical speed of a theoretical fracture filling type hydrate layer by utilizing a porosity data body and the anisotropic rock physical model;
and S4, comparing the inversion speed with the theoretical calculation speed to estimate the saturation of the fracture filling type hydrate.
2. The method of anisotropic seismic saturation prediction of a fracture-filled natural gas hydrate layer according to claim 1, wherein the synthetic seismic record comprises:
calculating a reflection coefficient sequence r (i) of different depths of the stratum by using stratum depth domain data of actual logging speed and density:
Figure FDA0004062943450000011
wherein the parameter is depth threshold, r (i) is the reflection coefficient of the ith layer, V i+1 、V i Logging speeds, p, for the i +1 th and i-th layers, respectively i+1 、ρ i Logging densities of an i +1 th layer and an i-th layer respectively;
and convolution is carried out with the input Rake wavelet information to generate a time domain synthetic seismic record:
Figure FDA0004062943450000012
in the formula, S (t) is synthetic seismic record of a time domain, d tau is a time interval, w (tau) is wavelets in offset trace concentration, and r (t-tau) is a time domain reflection coefficient sequence after velocity conversion of a depth domain reflection coefficient sequence r (i).
3. The method of anisotropic seismic saturation prediction of fracture-filled natural gas hydrate layers according to claim 1, wherein the post-stack seismic wave impedance inversion comprises: well logging data processing, seismic data processing, synthetic seismic recording, well seismic calibration, low-frequency model establishment, model analysis and updating.
4. The method of predicting anisotropic seismic saturation of a fracture-filled natural gas hydrate layer according to claim 1, wherein the calculating of the velocity, density and porosity data volume comprises: according to the velocity and density low-frequency model obtained by the inversion of the post-stack seismic wave impedance, dividing the velocity and density low-frequency model by the wave impedance data volume to obtain a density and velocity data volume; and calculating a porosity data volume using the density data volume:
Figure FDA0004062943450000021
where phi is the porosity, rho s Is the particle density, ρ is the density obtained by inversion, ρ w Is the density of water.
5. The method of anisotropic seismic saturation prediction of a fracture-filled natural gas hydrate layer according to claim 1, wherein the step S3 comprises:
calculating a Lame coefficient of a saturated water stratum by using a particle contact model by taking a porosity data volume and regional stratum lithology as input parameters;
and calculating the speed of the fracture filling type hydrate layer by using Backus average by taking the Lame coefficients of the saturated water stratum and the hydrate as input parameters.
6. The method of anisotropic seismic saturation prediction of a fracture-filled natural gas hydrate layer according to claim 5, wherein the particle contact model is used to calculate saturated water formation elastic parameters, and the formula comprises:
bulk modulus K of dry skeleton of sediment dry And shear modulus G dry ，
Figure FDA0004062943450000022
Or
Figure FDA0004062943450000031
Wherein φ is porosity; phi is a c Is the porosity of the dense random inclusion, about 40%; k HM And G HM Needs to be obtained by Hertz-Mindlin theoretical calculation; calculating the bulk modulus K and the shear modulus G of the solid particles according to mineral components by using a Hill average formula;
bulk modulus K of saturated water deposit sat And shear modulus G sat Calculated using the Gassmann equation:
Figure FDA0004062943450000032
G sat ＝G dry ；
in the formula, K w Is the bulk modulus of water;
calculating the Lame coefficient lambda and mu of the saturated water deposit according to the formula:
λ＝K sat -2μ/3,
μ＝G sat .
7. the method for predicting anisotropic seismic saturation of a fracture-filled natural gas hydrate layer according to claim 5, wherein the Backus is divided equally into two end members of a saturation stratum and a hydrate, and the calculation formula comprises:
Figure FDA0004062943450000033
Figure FDA0004062943450000034
Figure FDA0004062943450000035
Figure FDA0004062943450000036
N＝<μ>,
ρ＝<ρ>,
Figure FDA0004062943450000037
in the formula, A, C, F, L, N, rho and Q are elastic parameters in a rigidity matrix respectively, lambda and mu are Lame coefficients of a saturated water formation and a hydrate respectively, and gamma is an included angle between a wave vector and a symmetry axis; the tip brackets represent the weighted average expression form:
<G>≡(η 1 G 1 +η 2 G 2 )
Figure FDA0004062943450000041
in the formula eta 1 And η 2 G corresponds to the calculation formula in the tip bracket in the above stiffness matrix formula for the volumes of the two end members;
calculating the longitudinal and transverse wave speeds of the crack filling type hydrate layer by using the elastic parameters,
Figure FDA0004062943450000042
Figure FDA0004062943450000043
Figure FDA0004062943450000044
in the formula, V p 、
Figure FDA0004062943450000045
And &>
Figure FDA0004062943450000046
The longitudinal wave velocity, the vertical polarization transverse wave and the horizontal polarization transverse wave are respectively, namely the theoretical calculation velocity.
8. The method of anisotropic seismic saturation prediction of a fracture-filled natural gas hydrate layer according to claim 1, wherein estimating the saturation of a fracture-filled hydrate comprises:
and calculating the error between the inversion speed and the theoretical calculation speed under different hydrate saturation conditions, wherein when the error is less than 1e-2, the saturation is the finally calculated saturation of the fracture filling type hydrate layer.
9. The anisotropic seismic saturation prediction device of the fracture filling type natural gas hydrate layer is characterized in that the anisotropic seismic saturation of the fracture filling type natural gas hydrate layer is predicted and estimated through the following program modules:
the seismic data synthesis module is used for generating seismic data of the synthetic seismic record by utilizing the velocity and density logging data through a convolution method;
the inversion and porosity data module is used for obtaining a seismic wave impedance data volume, a velocity and density low-frequency model of the fracture filling type hydrate layer by using seismic data, velocity and density logging data and performing post-stack seismic wave impedance inversion, and calculating the velocity, density and porosity data volume of the hydrate layer;
the theoretical velocity calculation module is used for constructing an anisotropic rock physical model of the fracture filling type hydrate by utilizing the Backus average and particle contact model and calculating the theoretical velocity of a theoretical fracture filling type hydrate layer by utilizing the porosity data volume and the anisotropic rock physical model;
the saturation prediction estimation module is used for comparing the inversion speed with the theoretical calculation speed to estimate the saturation of the fracture filling type hydrate;
and the visualized mapping module is used for mapping each data volume in the process of predicting and estimating the saturation degree and visually displaying the data volumes.
10. Anisotropic seismic saturation prediction apparatus for a fracture-filled natural gas hydrate layer, comprising a processing unit and a storage unit, wherein the storage unit stores a program, and the processing unit loads the program to perform the method steps according to any one of claims 1 to 8 to predict and estimate the anisotropic seismic saturation of the fracture-filled natural gas hydrate layer.
CN202210398462.8A 2022-04-15 2022-04-15 Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer Active CN114879259B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210398462.8A CN114879259B (en) 2022-04-15 2022-04-15 Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210398462.8A CN114879259B (en) 2022-04-15 2022-04-15 Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer

Publications (2)

Publication Number Publication Date
CN114879259A CN114879259A (en) 2022-08-09
CN114879259B true CN114879259B (en) 2023-03-31

Family

ID=82670434

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210398462.8A Active CN114879259B (en) 2022-04-15 2022-04-15 Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer

Country Status (1)

Country Link
CN (1) CN114879259B (en)

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9482771B2 (en) * 2012-11-05 2016-11-01 Fugro Marine Geoservices, Inc. Method of indicating the presence of gas hydrate and shallow gas in deepwater environment
US9733372B2 (en) * 2012-11-23 2017-08-15 Fugro Marine Geoservices, Inc. Method and system for identification of gas hydrates and free gas in geologic beds
CN111722282B (en) * 2020-06-18 2022-11-29 中国科学院海洋研究所 Method for predicting natural gas hydrate reservoir top hydrate saturation by AVO
CN112149282B (en) * 2020-08-28 2024-05-28 中国石油天然气集团有限公司 Rock physical calculation method and system for saturation of natural gas hydrate in well
CN113376709B (en) * 2021-06-21 2022-05-13 西南石油大学 Method for predicting reservoir natural gas hydrate saturation by using logging data
CN113466963B (en) * 2021-08-13 2022-02-01 广州海洋地质调查局 Fracture filling II type hydrate saturation estimation method and processing terminal

Also Published As

Publication number Publication date
CN114879259A (en) 2022-08-09

Similar Documents

Publication Publication Date Title
CA3122509C (en) Machine learning-augmented geophysical inversion
US10572611B2 (en) Method and system for characterizing fractures in a subsurface region
US6993433B2 (en) Modeling gravity and tensor gravity data using poisson&#39;s equation for airborne, surface and borehole applications
US8352190B2 (en) Method for analyzing multiple geophysical data sets
AU2005250760A1 (en) Method for predicting lithology and porosity from seismic reflection data
CN104950331A (en) Earthquake prediction method for porosity and shale content of sand shale reservoir
CA2828425A1 (en) Sensitivity kernel-based migration velocity analysis in 3d anisotropic media
CN103513277B (en) Seismic stratum fracture crack density inversion method and system
CN105388520B (en) Seismic data prestack reverse time migration imaging method
CN103293553B (en) Upper and lower cable seismic acquisition data boundary element continuation bearing calibration at the bottom of a kind of Complex Sea
AU2019237361B2 (en) System and method for assessing the presence of hydrocarbons in a subterranean reservoir based on seismic inversions
US20190129049A1 (en) System and method for assessing the presence of hydrocarbons in a subterranean reservoir based on time-lapse seismic data
CN111505714B (en) Elastic wave direct envelope inversion method based on rock physical constraint
CN103907032A (en) Seismic imaging systems and methods employing correlation-based stacking
Li et al. A rock-physical modeling method for carbonate reservoirs at seismic scale
CN112147677B (en) Method and device for generating parameter tag data of oil and gas reservoir
CN114910950B (en) Method and system for estimating saturation of free gas below natural gas hydrate reservoir
CN115857006B (en) Submarine acoustic and physical parameter detection method, medium and system
CN114879259B (en) Method for predicting anisotropic seismic saturation of fracture filling type natural gas hydrate layer
CN109738944B (en) Wide-angle reflection-based seismic acquisition parameter determination method and device
CN115586572B (en) Seismic rock physical analysis inversion method for pore parameters and reservoir parameters
WO2012060888A1 (en) System and method for providing a physical property model
CN111175822B (en) Strong scattering medium inversion method for improving direct envelope inversion and disturbance decomposition
CN104483702A (en) Earthquake forward modeling method applicable to non-uniform motion water bodies
AlKawai et al. Integrating statistical rock physics and pressure and thermal history modeling to map reservoir lithofacies in the deepwater Gulf of Mexico

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
CP03 Change of name, title or address
CP03 Change of name, title or address

Address after: No. 7 Nanhai Road, Shinan District, Qingdao City, Shandong Province

Patentee after: Institute of Oceanology Chinese Academy of Sciences

Country or region after: China

Patentee after: OCEAN University OF CHINA

Patentee after: Qingdao Marine Science and Technology Center

Address before: 266071 Shandong province Qingdao City Nanhai Road No. 7

Patentee before: Institute of Oceanology Chinese Academy of Sciences

Country or region before: China

Patentee before: OCEAN University OF CHINA

Patentee before: QINGDAO NATIONAL LABORATORY FOR MARINE SCIENCE AND TECHNOLOGY DEVELOPMENT CENTER