Table Of ContentPublications of the University of Miskolc, Series A, Mining, Vol. 59. Geosciences, (2001), pp.115-138
THE INVERSION OF WELL LOG DATA USING SIMULATED ANNEALING
METHOD
M. Dobróka1-2 and P. N. Szabó2
Abstract
In the paper we present a global optimization method called Simulated An-
nealing is used to solve the nonlinear geophysical well logging inverse problem.
There are two possible ways to solve this problem. The first one is a conventional
inversion method, which estimates the unknowns in different depth points sepa-
rately. The other way is to get these parameters is the so called interval inversion
procedure, which uses the data set along the whole interval in a joint inversion pro-
cess. We made numerical investigations by synthetic and in-situ well log data to
compare these two algorithms so as to decide which is more powerful and yields
more accurate information about the downhole geological situation.
1. Introduction
The first stage of petrophysical interpretation for hydrocarbon exploration
is the evaluation of data sets measured in bore-holes. This operational problem is
generally solved by modern inversion methods with the application of up-to-date
informatics. The main object of this procedure is the detection of reservoir rocks in
the succession of strata, and determination of their layer-thicknesses and some pet-
rophysical parameters. Consequently, the unknowns of the inverse problem are
1 MTA-ME Research Group of Geophysical Inversion and Tomography 3515 Miskolc-Egyetemváros
2 Department of Geophysics, University of Miskolc, 3515 Miskolc-Egyetemváros
115
M. Dobróka and P. N. Szabó
geometrical parameters and several petrophysical characteristic values with rela-
tively lot of zoneparameters of the formations, where the petrophysical parameters
have constant values in a layer, but the zoneparameters are invariable in a longer
depth-interval. These physical and geometrical parameters are also used in hydro-
carbon forecasting and exploitation planning, therefore it is essential to determine
them as accurate as possible. In the industrial practice there are several softwares
for this goal, e.g. the CLASS, SAND, CRA which are deterministic procedures, but
the OPTIMA, OPTIMA 10 are probabilistic methods making complex reservoir
analyses by processing field data based on linear optimization methods. It is well-
known, that the linearized inversion methods are relatively quick, but the procedure
can often be trapped in local minima. Because of this fact the use of global optimi-
zation methods is reasonable.
2. The forward problem
In formulating the forward problem let us introduce the transposed vector
of the petrophysical model parameters as
m ={ POR,SXO,SW ,VCL,VSD,VLM }T (1)
where
POR - porosity,
SXO - water saturation in flushed zone,
SW - water saturation beyond invaded zone,
VCL - specific volume of clay,
VSD - specific volume of sandstone,
VLM-specific volume of limestone.
116
The inversion of geophysical well logging data...
To determine these parameters we utilise logs measured in the drill-hole, that rec-
ord different parameters of natural and/or induced physical fields in the function of
depth, yielding continuos geological information on the formations of the drilled
hole. Let us take the following transposed vector, which contains the observed data
of a possible combination of well logs:
d(oil) = {SP, GR, PORN, DEN, A T, RMLL, RLLD}7 (2)
where
SP - spontaneous potential data,
GR - natural gamma-ray data,
PORN- neutron porosity data,
DEN- density (scattered gamma-ray) data,
AT- acoustic traveltime data,
RMLL - microresistivity data,
RLLD - macroresistivity data.
The measured data reflect the immediate vicinity of the downhole, the SP and GR
logs are mainly sensitive to Iithology, the PORN-DEN-AT logs indicate porosity,
and RMLL-RLLD data are influenced by water saturation. Since the measurements
are carried out in a relatively complicated borehole surroundings both in physical
and geometrical sense (high pressure and temperature, vertically and horizontally
inhomogeneous layers, anisotropy, changing hole calliper, different penetration-
depth and vertical resolution capability of tools etc.), therefore at first we need to
correct these data then match these logs in pursuance of depth. In the first step of
forward modelling we create a petrophysical model for the formation of interest on
the basis of the observed data set and a priori information (e.g. lab-examination of
cores), characterising quantitatively each layer. In the next step we calculate physi-
cal data from the parameters of this model by certain petrophysical relationships.
117
M. Dobróka and P. N. Szabó
These relations are called response functions and they associate the model pa-
rameters with log data. The direct problem is solved by these multivariable func-
tions, where there is a g connection between the predicted model parameter vector
and the approximated data vector by calculation
dtmk)=g(m,c), (3)
where c- vector of set values like DENSD, ATCL etc. Let us survey the detailed
response functions of (2), which is applied also in a well logging data interpretation
system named ULTRA by Halliburton-Gearhart Co. (TH means theoretical value):
SPTH = SPSD - VCL (SPSD - SPSH), (4)
DENTH = POR[DENMF -C(l-SXO)] + VCL DENCL
" (5)
+ £ VMA DEN MA,
j
i=l
GRTH - {POR[GRMF SXO DENMF + GRCH{l - SX0)DENCH]
+ VCL • GRCL • DENCL + £ VMA, DENMA, • GRMA,} (6)
i=i
/ DENTH
PORNTH = POR[PORNMF - BCOR( 1 - SXO ) + BC] +
n (7)
VCL • PORNCL + £ VMA • PORNMA,
i
ATTH = PORfSXO ATMF + ATCH (1-SXO)] + VCL ATCL
(8)
+ YVMA ATMA
j i jt
f nN_~\
VCL ' POR 2 i
+ SXO 2 , (9)
yjRCL siBA RMF
/ r x o th
118
The inversion of geophysical well logging data...
1
(10)
•JRTTH
where n - the number of mineral components of rocks. On the left side of Eq. (4)-
(10) we can find the calculated data in bold letter, and on the right the model pa-
rameters stand in the same type, and there are also a lot of constants. Obviously,
the set of equations is nonlinear referring to the model parameters, so it is better to
solve the inverse problem by the help of a global optimalization method. During
the forward modelling we substitute the initial (and later the estimated) values of
(1) model parameters into the Eq. (4)-(10) equations, then after having data re-
ceived by this procedure we compare them with in-situ data and make a prediction
for the real petrophysical model by an inversion method.
3. Inversion algorithms
The nonlinear well logging inverse problem is conventionally solved by
separated inversion which estimates the petrophysical model parameters in a depth-
point by the use of all available data in the point measured by different physical
principles. We allow for that adjacent depth-points independent, this way we can
not determine the layer thicknesses and we must fix the zoneparameters during the
procedure as the inverse problem is marginally overdetermined. The theoretical
data in the point is calculated by means of a local set of response equations (see Eq.
(4)-(10)). The calculated data vector of the y-th log can be written in a general
forms as
d^=gJ{ml,...,mu), (11)
where M- the number of model parameters in the depth-point.
119
M. Dobróka and P. N. Szabó
Since the number of observed data is only just more than the unknown
model parameters in the point, this narrow overdetermination causes that the accu-
racy and reliability of the estimation of model parameters is relatively limited. For
example there are six petrophysical parameters in (1) against seven well log data in
(2). Therefore it is worth inverting data of a greater interval jointly in one proce-
dure. The so called interval inversion algorithm is based on the series expansion of
the petrophysical parameters. It develops depth dependent layer characteristic pa-
rameters DOBRÓKA (1991):
(12)
where m - the z'-th model parameter, z - the depth, B^1 - series expansion coeffi-
cients, - basis functions. For instance for layerwise homogeneous model we
can describe the model with the combination of unit step functions adopting the
least unknown parameters:
0, if z < z_
q
1 , if z„ . < z < z
' q -1 i
0 , if z > z
(If the petrophysical parameter shows vertical variation in a layer, we can choose
power or other appropriate basis function types). With this development in a series
the (11) connection modifies, it is now a response function interpreted in a depth-
interval. The synthetic data calculated from they'-th log is
120
The inversion of geophysical well logging data...
(13)
where the unknown parameters are the B coefficients. By optimalization we can
determine B-s that approach the petrophysical model parameters along the overall
observed interval after substituting them into (12) formula.
3.1. Linearized inversion methods
These conventional search methods reduce the solution of the nonlinear in-
verse problem to sequences of linear problems. Before using these procedures, the
problem should be linearized. Let m be a point in the parameter space close to the
0
solution of the nonlinear inverse problem
m-m +őm , (14)
0
where 8 m denotes the parameter correction vector. Solving the linearized inverse
problem means that in a method of successive approximation m is chosen m in
0
the subsequent iteration step, and the process will be continued until a stop crite-
rion is satisfied.
3.1.1. Linearized separated (point by point) inversion
Let us assume that the M-dimensional model parameter column vector and
the Z-dimensioned transposed calculated data vector are:
in ={m,m ,...,m } T
l 2 M
d ={d,d ,...,d }T.
t 2 L
121
M. Dobróka and P. N. Szabó
For an overdetermined inverse problem our goal is to determine the elements of the
model parameter vector by this equation
d = g(m), (15)
where the functional dependence is generally nonlinear, and the number of the data
is greater then those of model parameters. Let us develop Eq. (15) in Taylor-series
in the neighbourhood of m :
0
d = dl0) +GSm ,
where G - Jacobi matrix. This lead us to a linear set of equations
Sd = GSm , (16)
where the normalized terms owing to the different order of magnitude of the well
log data are
m.
jfobn)
8d =- a"
k (obs )
d
k
p m (0) 'dd<k>^
V DM, Y
Whereas Eq. (16) is inconsistent and the existence of measurement and model er-
rors the prediction error vector of measured and calculated data is different from
122
The inversion of geophysical well logging data...
zero in any case. Therefore the problem can be solved by the minimalization by the
following norm
(17)
which makes the parameter correction vector in every iteration step during the (It-
eratively Reweighted Least Squares) optimization. This is in the 17-th step:
őm^ =(^rW(9)GyGTW('l)Sd
p-2
M
W =
The goodness of inversion results are characterized in every depth-points sepa-
rately with the relative model and data distances. They are in case ofp = 2 (Least
Squares Method)
1 L (W'"*1) /ft"*)
100%
1(obx)
V
V M ) - mi'""0
Dmod = J— Y •100% (18)
(exact)
'Mtr m;
where m(esl) - the estimated model vector, and mlaaa) - the exact model in case of
inversion of synthetic data sets. The estimation error of the petrophysical parame-
ters and the reliability of the solution can be characterized by the elements of the
covariance and correlation matrixes (MENKE, 1984).
123
M. Dobróka and P. N. Szabó
3.1.1. Linearized interval inversion
Like in 3.1.1. the depth dependent layer parameters can be determined by means
of the minimalization the p-norm of the error vector of the overall observed and
calculated well log data:
DP L I
• I
Sd* (19)
-
/A
where DP - the number of depth-points in the whole observed interval. The good-
ness of result are
DP l f j _/y('"/c) V
= J -1 J —DP ftL 14 'hhi i w hi 100%
Ddata »(oA.v)
VDPLtijÍ
V >" /
5
•100%, (20)
{exact)
rk
where S- the number of the layers in the interval.
3.2. Global optimization method - Simulated Annealing
The 3.1. methods are the most prevailing in the practice of inversion, because
in the case of having an initial model nearby the solution they are very quick and
effective algorithms. They are also capable of quality checking of the estimated
model parameters. After all they aren't absolute minimum searching methods, thus
they assign the solution to a local optimum of the (17) or (19) objective function
124
Description:compare these two algorithms so as to decide which is more powerful and yields and RMLL-RLLD data are influenced by water saturation.
