Method for simultaneous inversion of five reservoir parameters of submarine gas reservoirs using P-WAVE velocity and density data
A method for simultaneous inversion of five reservoir parameters of submarine gas reservoirs using P-wave velocity and density data, comprising: acquiring an observed P-wave velocity, an observed density, and a reservoir distribution range and constructing a P-wave velocity model and a bulk density model. Obtaining variable ranges of five reservoir parameters from geological experience, and setting unknown variables to be inverted, and constructing an optimization problem. Using an interior-point algorithm to iteratively solve the optimization problem, outputting optimal inversion values and consolidation factors for the five reservoir parameters. Obtaining the 2D/3D P-wave velocity profile and a 2D/3D density profile; using the optimal consolidation factor, inputting them into the optimization problem constructed with the P-wave velocity model and the bulk density model, repeating the inversion process for each CDP trace to invert the five reservoir parameters for each depth point.
1 . A computer-implemented method for simultaneous inversion of five reservoir parameters of submarine gas reservoirs using P-wave velocity and density data, comprising:
acquiring, via a logging equipment, an observed P-wave velocity
(
V
p
obs
)
,
an observed density
(
ρ
b
obs
)
,
and a reservoir distribution range, and measuring an actual value of a free gas saturation, a porosity, and mineral composition ratios of a rock skeleton to validate inversion accuracy; constructing a physically-based P-wave velocity model and a bulk density model of the a reservoir rock using a simplified two-phase Biot equation that incorporates a depth-dependent consolidation factor (α) which is depth-dependent, a weighted average equation of mixed fluid that uses a calibration factor (χ) to model patchy saturation, and a volume density equation, wherein these models incorporate physical properties of clay, quartz, carbonate, sea water, and natural gas;
obtaining experienced range values of a free gas saturation, a range of the porosity, and a range of each of the mineral composition ratios; setting the free gas saturation (S g ), the porosity (φ), a clay content (V clay ), a carbonate content (V car ), and a water-gas mixture calibration factor (χ) as five unknown variables to be inverted in the P-wave velocity model and the bulk density model; associating the P-wave velocity model with the observed P-wave velocity and the bulk density model with the observed density; and
constructing a multi-dimensional nonlinear constrained optimization problem with a goal of minimizing a sum of mean squared errors between model-calculated values and observed values, while setting physically-based variable ranges for each inverted unknown wherein, 0≤S g ≤0.5, 0≤V clay ≤1,
0≤ V car ≤1, 0≤ V clay +V car ≤1, 0≤χ≤1;
after constructing the optimization problem, simultaneously obtaining values of the five reservoir parameters from logging data and initial iteration values, by using multiple sets of initial iteration values randomly selected within the variable ranges and solving the optimization problem with an interior-point algorithm, while dynamically adjusting the consolidation factor for each of the reservoir distribution points during iteration a solution process; and judging inversion accuracy by calculating a root mean square error between inversion results and the observed reservoir parameters, until convergence; and outputting inversion values for validating inversion accuracy together with a corresponding consolidation factor used for predicting a two dimensional (2D) reservoir parameter prediction or a three dimensional (3D) reservoir parameter;
during a prediction process, acquiring submarine gas reservoir seismic data collected by a seismic equipment, inverting the submarine gas reservoir seismic data, using a post-stack inversion method for post-stack seismic data or an amplitude variation with offset (AVO) inversion method for pre-stack seismic data, to obtain a 2D P-wave velocity profile or a 3D P-wave velocity profile and a 2D density profile or a 3D density profile;
using the corresponding consolidation factor, inputting the 2D P-wave velocity profile or 3D P-wave velocity profile, the 2D density profile or the 3D density profile as inputs to the optimization problem constructed with the P-wave velocity model and the bulk density model, repeating an inversion process of each seismic common depth point (CDP) trace and each depth point to obtain the inverted unknown variables of the five reservoir parameters; outputting a spatially indexed 2D reservoir parameter distribution profile or a spatially indexed 3D reservoir parameter distribution profile of the submarine gas reservoir, thereby enabling improved reservoir volume evaluation and development optimization through spatially-distributed physical parameters of the gas reservoir.
2 . The method of claim 1 , wherein the bulk modulus K s of a rock skeleton mineral matrix and the shear modulus us of the rock skeleton mineral matrix are calculated by a Hill formula, and the Hill formula is as follows:
K
S
=
{
[
V
clay
K
clay
+
(
1
-
V
clay
-
V
car
)
K
q
u
a
+
V
car
K
c
a
r
]
+
1
[
V
clay
/
K
clay
+
(
1
-
V
clay
-
V
car
)
/
K
q
u
a
+
V
car
/
K
car
]
}
/
2
;
and
μ
S
=
{
[
V
clay
μ
clay
+
(
1
-
V
clay
-
V
car
)
μ
q
u
a
+
V
c
a
r
μ
c
a
r
]
+
1
[
V
clay
/
μ
clay
+
(
1
-
V
clay
-
V
car
)
/
μ
q
u
a
+
V
c
a
r
/
μ
c
a
r
]
}
/
2
;
wherein, V clay is a volume fraction of the clay in the rock skeleton mineral matrix, and V car is a volume fraction of the carbonate in the rock skeleton mineral matrix.
3 . The method of claim 2 , wherein the bulk modulus approximate Biot coefficient β P and shear modulus approximate Biot coefficient β S are as follows:
β
P
=
ϕ
(
1
+
α
)
(
1
+
αϕ
)
;
β
S
=
ϕ
(
1
+
γα
)
(
1
+
γαϕ
)
;
wherein
γ
=
(
1
+
2
α
)
(
1
+
α
)
,
α
=
α
0
(
d
0
d
)
1
/
3
;
α is a consolidation factor at a target point, α 0 is a consolidation factor at a reservoir bottom, d 0 is a depth of reservoir bottom, d is a target depth, and φ is the porosity.
4 . The method of claim 3 , wherein a variable range of the consolidation factor at the reservoir distribution point is 1≤α 0 ≤15.
5 . The method of claim 1 , wherein the mineral composition ratios comprise a clay bulk modulus K clay , and a clay shear modulus μ clay , and a quartz bulk modulus K qua , and a quartz shear modulus μ qua , and a carbonate bulk modulus K car , and a carbonate shear modulus μ car .
6 . The method of claim 1 , wherein the simplified two-phase Biot equation is as follows:
K
=
K
S
(
1
-
β
P
)
+
β
P
2
K
av
,
K
av
=
1
/
(
β
P
-
ϕ
K
S
+
ϕ
K
f
)
,
μ
=
μ
S
(
1
-
β
S
)
;
wherein, K is an overall bulk modulus of a rock, K s is a bulk modulus of a rock skeleton mineral matrix, μ S is a shear modulus of the rock skeleton mineral matrix, β P is a bulk modulus approximate Biot coefficient, K av is an equivalent bulk modulus, φ is the porosity, K f is a water-gas mixed fluid bulk modulus, μ is the overall shear modulus of the rock, μ S is a shear modulus of the rock skeleton mineral matrix, and β S is a shear modulus approximate Biot coefficient.
7 . The method of claim 1 , wherein, according to the weighted average equation of mixed fluid, the water-gas mixed fluid bulk modulus K f is expressed as follows:
K
f
=
χ
[
(
1
-
S
g
)
K
w
+
S
g
K
g
]
+
(
1
-
χ
)
/
[
(
1
-
S
g
)
/
K
w
+
S
g
/
K
g
]
;
wherein χ is a water-gas mixture calibration factor, and 0≤χ≤1, S g is the free gas saturation, K w is a bulk modulus of liquid water, and K g is a bulk modulus of free gas.
8 . The method of claim 1 , wherein the bulk volume equation is as follows:
ρ
b
=
(
1
-
ϕ
)
[
V
clay
ρ
clay
+
(
1
-
V
clay
-
V
car
)
ρ
qar
+
V
car
ρ
car
]
+
(
1
-
S
g
)
ϕρ
w
+
S
g
ϕρ
g
;
wherein, ρ b is a reservoir bulk density, V clay is a volume fraction of the clay, ρ clay is a density of the clay, ρ qar is a density of the quartz, ρ car is a density of the carbonate, ρ w is a density of the liquid water, and ρ g is a density of the free gas.
9 . The method of claim 8 , wherein the reservoir P-wave velocity in the P-wave velocity model is as follows:
V
P
=
K
+
4
μ
/
3
ρ
b
;
wherein, V P is a reservoir P-wave velocity, K is the overall bulk modulus of a rock, and μ is an overall bulk modulus of the rock matrix.
10 . The method of claim 1 , wherein an objective function corresponding to multi-dimensional nonlinear constrained optimization problem is as follows:
min
f
(
S
g
,
ϕ
,
V
c
l
a
y
,
V
c
a
r
,
χ
)
=
[
V
p
o
b
s
-
V
p
m
o
d
(
S
g
,
ϕ
,
V
c
l
a
y
,
V
c
a
r
,
χ
)
]
2
+
[
ρ
b
o
b
s
-
ρ
b
m
o
d
(
S
g
,
ϕ
,
V
c
l
a
y
,
V
c
a
r
)
]
2
;
wherein, f is the objective function,
V
p
o
b
s
is the observed P-wave velocity,
V
p
m
o
d
is a calculated P-wave velocity from the P-wave velocity model,
ρ
b
o
b
s
is the observed density, and
ρ
b
m
o
d
is a calculated bulk density from the bulk density model.
11 . The method of claim 1 , wherein the five reservoir parameters comprise the free gas saturation, the porosity, a clay volume fraction, and a carbonate volume fraction from the logging data.
12 . The method of claim 1 , wherein the root mean square error comprises a first root mean square error between inverted free-gas saturation values and corresponding well-log free-gas saturation values and a second root mean square error between inverted porosity values and corresponding well-log porosity values.
13 . The method of claim 1 , wherein when solving the optimization problem using the interior-point algorithm, the multiple sets of initial iteration values are randomly selected within the variable range, and the result with the smallest objective function value is selected as the corresponding consolidation factor after solving each set.
14 . The method of claim 1 , wherein the step of repeating the inversion process comprises:
repeating the following steps: constructing the P-wave velocity model and the bulk density model based on mineral physical property data; setting the free gas saturation, the porosity, the clay content, the carbonate content, and the water-gas mixture calibration factor as unknown variables to be inverted, and constructing the optimization problem with a goal of minimizing the sum of squares of the model calculation values and profile data; using the interior-point algorithm for iterative solving and outputting the five reservoir parameters for each seismic CDP trace depth point.
15 . The method of claim 1 , wherein the 2D reservoir parameter distribution profile or the 3D reservoir parameter distribution profile comprises: spatial distribution information of the free gas saturation, the porosity, the clay content, the carbonate content, and the water-gas mixture calibration factor in the submarine gas reservoir region; the spatial distribution information corresponds one to one with each seismic CDP trace and each depth point.
16 . A system for simultaneous inversion of five reservoir parameters of a submarine gas reservoir using P-wave velocity and density data, comprising a memory, a processor, and a program for simultaneous inversion of a five reservoir parameters of submarine gas reservoirs using P-wave velocity and density data, wherein the program is executed by the processor to perform steps of the method for simultaneous inversion of five reservoir parameters of the submarine gas reservoirs based on the P-wave velocity and the density of claim 1 .
17 . A non-transitory computer-readable storage medium, wherein the non-transitory computer-readable storage medium is configured for storing a program for simultaneous inversion of five reservoir parameters of submarine gas reservoirs using P-wave velocity and density data; when the program is executed by a processor to perform steps of the method for simultaneous inversion of five reservoir parameters of the submarine gas reservoirs based on the P-wave velocity and the density of claim 1 .
18 . A device for simultaneous inversion of five reservoir parameters of submarine gas reservoirs using P-wave velocity and density data, comprising:
data acquiring and rock physics modeling module, acquiring, via a logging equipment, an observed P-wave velocity
(
V
p
obs
)
,
an observed density
(
ρ
b
obs
)
,
and a reservoir distribution range, and measuring an actual value of a free gas saturation, a porosity, and mineral composition ratios of a rock skeleton to validate inversion accuracy; constructing a physically-based P-wave velocity model and a bulk density model of the a reservoir rock using a simplified two-phase Biot equation that incorporates a depth-dependent consolidation factor (α) which is depth-dependent, a weighted average equation of mixed fluid that uses a calibration factor (χ) to model patchy saturation, and a volume density equation, wherein these models incorporate physical properties of clay, quartz, carbonate, sea water, and natural gas;
geological experience obtaining and optimization problem constructing module, obtaining experienced range values of a free gas saturation, a range of the porosity, and a range of each of the mineral composition ratios; setting the free gas saturation (S g ), the porosity (φ), a clay content (V clay ), a carbonate content (V car ), and a water-gas mixture calibration factor (χ) as five unknown variables to be inverted in the P-wave velocity model and the bulk density model; associating the P-wave velocity model with the observed P-wave velocity and the bulk density model with the observed density; and
constructing a multi-dimensional nonlinear constrained optimization problem with a goal of minimizing a sum of mean squared errors between model-calculated values and observed values, while setting physically-based variable ranges for each inverted unknown variable; wherein, 0≤S g ≤0.5, 0≤V clay ≤1,
0≤ V car ≤1, 0≤ V clay +V car ≤1, 0≤χ≤1;
optimization solution and parameter outputting module, after constructing the optimization problem, simultaneously obtaining values of the five reservoir parameters from logging data and initial iteration values by using multiple sets of initial iteration values randomly selected within the variable ranges and solving the optimization problem with an interior-point algorithm, while dynamically adjusting the consolidation factor for each of the reservoir distribution points during iteration; and judging inversion accuracy by calculating a root mean square error between inversion results and the observed reservoir parameters, until convergence; and
outputting inversion values for validating inversion accuracy together with a corresponding consolidation factor for predicting a two dimensional (2D) reservoir parameter prediction or a three dimensional (3D) reservoir parameter;
seismic data inversion and reservoir parameters distribution inversion outputting module, during a prediction process, acquiring submarine gas reservoir seismic data collected by a seismic equipment, inverting the submarine gas reservoir seismic data, using a post-stack inversion method for post-stack seismic data or an amplitude variation with offset (AVO) inversion method for pre-stack seismic data, to obtain a 2D P-wave velocity profile or a 3D P-wave velocity profile and a 2D density profile or a 3D density profile;
using the corresponding consolidation factor, the 2D P-wave velocity profile or 3D P-wave velocity profile, the 2D density profile or the 3D density profile as inputs to the optimization problem constructed with the P-wave velocity model and the bulk density model, repeating an inversion process of each seismic common depth point (CDP) trace and each depth point to obtain the inverted unknown variables of the five reservoir parameters; outputting a spatially indexed 2D reservoir parameter distribution profile or a spatially indexed 3D reservoir parameter distribution profile of the submarine gas reservoir, thereby enabling improved reservoir volume evaluation and development optimization through spatially-distributed physical parameters of the gas reservoir.