Reflection of elastic waves in the layered Biot medium

A.A. Vikhorev1 , Mike Ammerman2 & E.M. Chesnokov3

1 Institute of the Physics of the Earth, Russian Academy of Science, Moscow, Russia;
2 Devon Energy Inc. OK, USA;
3 Sarkeys Energy Center University of Oklahoma, OK, USA

ABSTRACT: A solution for the generalized system of Biot equations has been obtained. It describes vertical
propagation of a P-wave in the multi-layered two-phase sedimentary rock model with provision for viscous
dynamic interaction between the phases. Based on the analysis of this solution the differences in the shape and
spectra of the pulses (seismic signals) reflected from the top and bottom bounds of the distinguished layer which
are possible in principle have been forecasted. This layer is distinguished from the others by the high liquid
concentration and has a good permeability. Special quantitative characteristics related to different viscosity
mechanisms for Newton and Maxwell liquids have been theoretically discovered.
Based on the short time spectral analysis of the velocity signals qualitative characteristics of the pulses
reflected from the top and bottom bounds of the productive layer have been defined. Within the Biot model the
above characteristics may be considered as hypothetical criteria for discovery of the productive layer on the basis
of seismic sounding.

1 INTRODUCTION of vertical sounding and profiling. These problems

attain significance in oil, gas and underground water
There is a class of direct problems which remains non- prospecting, where it is required to determine the
elementary but still allows a possibility of analytical principle existence of the productive layer and make
study and may reveal qualitative characteristics of the some quantitative estimates. At the moment the
mediums structure. In particular they allow to reveal AVD (Absorption and Velocity Dispersion) method
the signs of presence or absence of a two-phase layer suggested by Rapoport, 19922000, Ryjkov, 1994 is
containing solid and liquid (gaseous) components, and successfully applied in oil and gas prospecting.
also signs of the existence of pores or cracks in a partic- The definition of the equation types which describe
ular layer within the multilayered medium. Therefore wave propagation in a medium saturated with liquid-
we consider a wave process in a multilayered model gas mixture in the state of phase equilibrium remains
where liquid saturation of pores and mean radius of a pressing issue. Model calculations also reveal the
pore channels, and thus dispersion and attenuation of frequency dependence of the reflection indices at the
waves could be set as unique for each layer. The size boundaries between the layers with different saturation
of the pore channels shows the significance of the and permeability.
effect of liquid motion on the overall wave motion This article is dedicated to a theoretical solution of
of the medium. It could be said that the radius of the wave sounding problem. Synthetic seismograms
pore channels and the viscosity of the filling liquid are obtained for a model which has an exact solution
determine the number of oscillational degrees of free- and defined by generalized Biot equations, includ-
dom involved in the wave process. Therefore the actual ing those accounting for the effects of the second
number of degrees of freedom may also vary from viscosity. The used method allows to calculate the
layer to layer. Non-standard dynamic properties of a seismograms and to observe the effects caused by
complex structured media, or in particular, frequency the presence of the distinguished layer for various
dependencies of the local properties form a foundation equation types and medium models. A comparison of
for possible detection of a distinguished layer within synthetic and experimental observations will allow to
a multilayered medium based on the wave sounding check the hypotheses of the mediums behavior and to
data. make the right choice of a qualitative theory. The latter
Observation and analysis of dispersion, attenua- will allow to extract more complete information from
tion and re-distribution of energy over the degrees of the real seismic measurements data, i.e. to solve the
freedom have more practical benefits in the process reverse problem within the chosen theory.

267
Copyright 2005 Taylor & Francis Group plc, London, UK
2 THE GENERALIZED WAVE EQUATION and the resultant pulse density through the compo-
FOR A TWO-PHASE MEDIUM nents (u0 , u). In the continuity equation (2) let us use
the linear relation between the adiabatic variations of
A two-phase medium consisting of an elastic perme- the density , the stress tensor ij and the porous liquid
able matrix (skeleton) and a viscous liquid is described pressure P. This relation follows from the general non-
by M. A. Biot equations (Biot, 1956) where the holonomic constitutive equation (Nikolayevsky, 1984)
dynamic interaction between the oscillating liquid and for this two-phase medium: F(, , P; , , S) = 0
P;
the matrix is set by the frequency function. This results where S is the entropy. Let us present the constitutive
in Biot medium showing the frequency dispersion equation in adiabatic variations and express the den-
of waves velocity and attenuation and possessing a sity variation. All of the variations shall be considered
higher number of oscillational degrees of freedom in as changes of the corresponding values in time, then
comparison with the standard elastic medium. The we will obtain the required linear dependence:
Biot system of equations (Biot, 1962) consists of
one scalar and two vector equations which define the
following conditions:
Motion equation of a small element of the medium: where

where is the resultant pulse density, is the resultant

stress tensor, P is the mean porous liquid pressure, f is and use it in the continuity equation (2). As a result
the outside force density. equations (1,2) form a system which is closed in rela-
Continuity equation: tion to the required functions {u (x, t); = 0, 1, 2, 3},

and introduction of a uniform denomination C (x; p)
for the coefficients of the obtained system of equations
will allow to put it in a universal form:
The resultant pulse density is the density of the mass
flow and its divergence is related to the velocity of
change in the overall mass density .
Motion equation of liquid relative to an element of where , , , = 0, 1, 2, 3 and partial time and coor-
the skeleton: dinates derivatives are written as 0 = /t, j = /xj .
Representation of the more complex mechanisms of
viscosity as processes with a limited relaxation time
leads to the appearance of the dynamic permeability
The relative liquid velocity wis proportionate to which shows resonant behavior at certain frequencies.
the gradient of pressure and density of the inertia force The latter is related to the transition from Newtonian
f . Viscous dynamic operator Y( p) describes the liquid to Maxwellian liquid which is characterized
frequency dependence of the viscous friction force. by the Deborah number . The corresponding type
Here u(x, t) is the displacement vector of the skeleton of viscous dynamic operator Y(p) has been obtained
points; w(x, t) = (uf u) is the mean displacement in the works by Tsiklauri, 2001 and further stud-
of the liquid relative to the skeleton, i.e. uf (x, t) is ied in the works of Zhi-Wen Cui, 2003. Supposing
the absolute liquid displacement; matrix porosity; that in the equation (3) p = i, where is the cir-
Cijkl is the resultant tensor of the elasticity mod- cular frequency let us introduce a denomination for
uli; = (1 )s + f is the macroscopic summary the transition frequency: c = /<a2 >, where is the
density made of densities of pure phases s , f . kinematical liquid viscosity and <a2 > is the mean
In order to put the system (13) in its universal form, squared radius of the porous channels. It should be
let us introduce a potential u0 (x, t) for the oscillating noted that in the static limit Y(0) = c is always true,
part of the pressure. In the system (13) we substitute: then equation (3) represents Darcys law. If Deborah
number has a finite value, then the variational deriva-
tive of density over the components of density
tensor ij with constant entropy S is a function of the
frequency:
where P0 is the equilibrium value of the mean porous
pressure. Then from the last equation (3) we can extract
the relative liquid velocity:

268
Copyright 2005 Taylor & Francis Group plc, London, UK
The general development of this result is given in existence/absence of a distinguished layer, for example
(Landau, Livshits, 1988). Here c0 and c are low the productive layer within a multi-layered medium.
and high frequency limits of the velocity parameter The universal equation (5) in this case may be
c. The characteristic process = 1/c is the relax- reduced to:
ation time of a certain process, for example a
phase transition or a process of energy re-distribution
between oscillational and rotational degrees of free-
dom in large molecules (heavy hydrocarbons, poly-
mers). Similar processes also called effects of sec- where  is a (4 4) matrix and where the wave field
ond viscosity occur in Maxwell liquids during the is a column:
wave passage. For Newtonian liquids process = 0, i.e.
the process does not occur and Deborah number
= (viscosity /process ) (approaches infinity). The
characteristic viscosity time viscosity = a2 / defines The components of the wave field here are as follows:
the transition frequency c = 1/viscosity in a two-phase 00 = 0 u0 = P/f , where P(x, t) is the oscillating
model of the heterogeneous medium. The total time part of the porous pressure; u0 (x, t) is the pressure
derivative contains partial derivatives over the coor- potential; 03 = 0 u3 is the velocity of skeleton points
dinates q = {, ij , P} which define the state of the oscillation, where u3 (x, t) is the vertical displacement
medium as a thermodynamic system. of skeleton points; 03 = C30 u , 3 = C3 u are the
3 3

components of the generalized stress tensor.

4 MODELS
Therefore the non-holonomic constitutive equation The simplest model of a horizontally homogeneous
describes a non-localized in time relation between Biot medium consists of two layers in half-space. All
stress and deformation which is characteristic of a the parameters are shown as graphs on Figure 1.
Maxwellian liquid and causes frequency dispersion. The models studied include one productive layer
This frequency dispersion is independent from the (PL) covered by one or more non-productive layers
effects of dissipation on inhomogenuities and from with low permeability. It is supposed that the effec-
the viscous dynamic interaction of solid and liquid tive transition frequency in the PL (collector rock)
phases. f0 = min{c , c }/2 is within the seismic frequency
range. The hypothesis that the value of the effective
transition frequency lies within the seismic range is
3 THE SOLUTION FOR A PLANE LAYERED confirmed by experimental studies of dispersion and
BIOT MEDIUM attenuation of P-waves based on vertical seismic pro-
filing (Ryjkov, 1994). Each model containing a PL is
The study of multi-layered model where dispersion set against a comparison model which differs from the
and attenuation of waves inside each layer are deter- former only by replacement of the productive layer
mined by the heterogeneous properties of the medium with a non-productive one.
presents a special interest. In particular the relative
volume of the liquid, its type, viscosity and the mean Depth, km
0.000
radius of porous channels may be set for each layer 0.100
Vs Vp Vp2 fo

along with the elastic properties of the medium and 0.200

75
0.300
layers thickness. The relative volume of the liquid 0.400
50

corresponds to the porosity of the matrix = (x) 0.500

30
15
while the mean radius of porous channels a = a(x) 0.600
0.700 2

accounts for the permeability of the matrix. Kinemat- 0.800

ical viscosity = (x) and Deborah number = (x)

0.900
1.000

are properties of the liquid. The transition frequency 1.100

1.200
defines the scale of the wave propagation frequency 1.300

dependence. In the studied model the above values are 1.400

0 10 20 30 2 3 4 0.0 0.5 1.0 0.0 0.3 10 100 100010000
set by the piecewise constant functions of the (verti-
K, ,
cal) coordinate x = x3 . Let us suppose that the push is GPa
Vp,Vs
GPa
Vp2
km/s
fo,
Hz
distributed along the upper bound (over the surface)
and generates a P-wave propagating in the vertical Figure 1. The simplest model of a plane layered two-phase
direction. The benefit of this model is that it pro- medium. Of the two layers in half-space the lower layer is
vides an exact solution which may reveal signs of productive i.e. saturated and permeable.

269
Copyright 2005 Taylor & Francis Group plc, London, UK
5 RESULTS arbitrary, but their thickness has a macroscopic value.
The spectral content of the pulses reflected from the
5.1 Signals of velocity and stress points top and bottom bounds shown on Figure 3 is the same
The solution obtained for the P-wave in a plane-layered (Fig. 4b).
medium simulates the seismic signals which in theory
could be observed through vertical sounding. Figure 2 5.3 The form of arrival pulses in principal cases
shows synthetic seismograms of velocity and stress in a
vertical P-wave calculated for two layers in half-space Further Figures 5, 6 and 7 show arrival pulses of the
(Fig. 1). P-wave reflected from the top and bottom bounds
These seismograms show the characteristic differ- of the distinguished layer. The signals V1(x, t) cor-
ences of the speed v = u 3 (x, t) and stress 33 (x, t) respond to displacement velocity in the productive
signals in the displacement pulses u3 (x, t) which are model. V1(0, t) is the velocity of vertical displace-
reflected from the top (1) and bottom (2) bounds of ment of surface points; V1(h, t) is the velocity of points
the PL. displacement at the depth h = 12 m. In the same way
The pulses are marked as 1 and 2 on Figure 2. These V2(x, t) are the velocity signals for the comparison
pulses are reflected correspondingly from the top and model. Figure 5 shows the form of reflected pulses in
bottom bounds of the PL. the case where the distinguished layer is saturated with
On the other hand Figure 3 shows the same seismo-
grams for the comparison model. Here the pulses are Model with productive layer
only different in amplitude and polarity which depends f0 = 5000 Hz
1 2
on the contrast of the layers.
(h,t)
.
u (h,t)
5.2 The short time spectral analysis
Figure 4 shows the results of short time spectral anal- u. (0,t)
ysis of the model with a permeable and saturated layer u. (h,t)
(Fig. 4a) and for the comparison model, where the layer u. (0,t) 
positioned at the same depth has low permeability and
u. (h,t)
low saturation. The said pulses 1 and 2 have differ-
ent spectral content (Fig. 4a). This result remains the 0.52 0.56 0.60 0.64
same even if the number of non-productive layers is Time, sec.

Figure 3. Comparison of signals on the surface and at the

signal of velocity on surface depth h = 12 m for a model with productive layer.
signal of velocity at depth h
signal of stress at depth h
Model with productive layer
f0 = 50 Hz
1 2
(h,t)
.
u (h,t)

.
u (0,t)
.
u (h,t)

.
u (0,t) 
.
u (h,t)

0.52 0.56 0.60 0.64

Time, sec. Figure 4. Short time spectral analysis of the wave sounding
signal: (a) The simplest model with one saturated and highly
Figure 2. Comparison of signals on the surface and at the permeable layer; (b) Corresponding result for the comparison
depth h = 12 m for a model with productive layer. (1 pulse model with the layer of low saturation and permeability. The
reflected from the top bound of the productive layer, 2 pulse red color shows the maximum level of spectral amplitude.
reflected from the bottom bound). The blue color is the minimum amplitude.

270
Copyright 2005 Taylor & Francis Group plc, London, UK
Newtonian liquid (e.g. water, Deborah number >> 1) liquid (Deborah number < 1) with the mean radius
with big radius of porous channels a = 50 mkm where of porous channels a = 10 mkm.
the transition frequency defined as fc = /a2 /2 is Here (Fig. 6) the transition frequency defined as
within the seismic frequency range. In this case its fc = /a2 /2 is outside of the seismic range 5000
value is 70 Hz. Then the spectral splitting of the first Hz, but the spectral splitting still occurs due to the
and the second pulses is observed. frequency dispersion caused by the processes occur-
The pulses on the Figure 6 are observed when ring in Maxwellian liquid during the wave passage.
the distinguished layer is saturated with Maxwellian The effective transition frequency when the Deborah
number < 1 is defined as c /2.
The models where the productive layer is simulta-
neously saturated with gas and liquid and the relative
gas volume is higher than the liquid represent a spe-
cial interest. In this case gas cannot be dissolved in
liquid and therefore there exists a phase boundary. The
scale of this boundary is comparable with the size of
porous channels. For example we can consider a local-
ized model of a porous channel coaxially filled with
gas and liquid. Naturally this will be stable if the liq-
uid is uniformly distributed over the channel walls. By
solving the border problem for a coaxial porous chan-
nel model we will obtain the type of function Y(p) for
the case of liquid and gas filling, where both gas and
liquid can have different values of Deborah number.
The result will be a three phase medium.
Figure 7 shows pulses reflected from the bounds of
such three phase medium. In the used model the radius
of porous channels is small a = 2 mkm, Deborah
number for the liquid is = 100 and for the gas

Figure 5. The productive layer saturated with Newtonian

liquid (big Deborah number) with big radius of porous
channels.

Figure 7. The productive layer is saturated with Maxwellian

gas (small Deborah number) and Newtonian liquid with
coaxial geometry of phase boundary in the porous channels
Figure 6. The productive layer saturated with Maxwellian (small porous channels radius).
liquid (small Deborah number) with medium radius of porous The general porosity is equal to 7.5%. This porous space is
channels. saturated with 75% gas and 25% liquid.

271
Copyright 2005 Taylor & Francis Group plc, London, UK
 Model with high-permeable layer boundary between the layers with different porous
saturated by Newton liquid & Maxwell gas channels radiuses or with different transition fre-
quencies. Apart from that an intensive absorption of
the high frequency components of the passing pulse
occurs in the saturated layer. The index of P-wave
refraction from the upper bound of the PL grows
with frequency, thus in the pulse 1 high frequency
harmonics are strengthened which leads to trans-
formation of the descending pulse and the number
of oscillations increases. On the other hand in the
(Arrival time (sec.)

pulse 2 which is reflected from the bottom bound

and has passed through the PL twice the high fre-
quency harmonics are strongly subdued. Therefore the
descending impulse also transforms. It becomes more
monotonous and lengthy, i.e. the number of oscilla-
tions decreases. The described spectral deviation of
the reflected pulses becomes more pronounced if Deb-
orah number becomes lower which has been studied
in the works of Tsiklauri, 2001, Zhi-Wen Cui, 2003.
Decrease of Deborah number characterizes transition
from Newtonian liquid to Maxwellian liquid.

7 CONCLUSIONS

1. An exact solution of the generalized system of Biot

equations has been obtained. It describes vertical
Figure 8. Short time spectral analysis of the wave sounding propagation of a P-wave in a horizontally layered
signal for a model with one three-phase layer saturated with two-phase medium.
gas and liquid. 2. Based on the short time spectral analysis of the
velocity signals the qualitative characteristics of
= 0.0001. Porosity equals 7.5%, 75% of which is pulses reflected from the top and bottom bounds of
taken by gas and 25% by liquid. the productive layer have been defined. The discov-
Here due to a relatively low gas density together ered spectral deviation of the reflected pulses can be
with the effects of second viscosity the absorption strictly defined as an algorithm of pulse selection.
spectrum is banded and there are frequency ranges Within the Biot model the above characteristics may
where attenuation of the second P-wave is low. Aside be considered as hypothetical criteria for productive
from that the frequency dependence of the wave prop- layer discovery on the basis of seismic sounding.
agation velocities is resonant. This leads to appearance
of an oscillating part of the second arrival (Fig. 7).
The resonant properties of the frequency depen-
dency of the wave velocity and attenuation lead to the REFERENCES
appearance of banded spectrum of the pulse which
went through the productive layer. In particular Fig- Biot M.A. 1956, Theory of propagation of elastic waves in
a fluid-saturated porous solid. J. Acoust.Soc.Am., vol. 28,
ure 8 shows the result of short time spectral analysis of p.168191.
the pulses presented in Figure 7. Therefore in the spec- Biot M.A. 1962, Mechanics of deformation and acoustic
trum of a broadband signal (wave) which went through propagation in porous media., J. Appl.Phys., vol. 33,
this reservoir we can expect to see the appearance p. 14821498.
of absorption lines and also propagation lines of cer- Dai N., Vafidis A., Kanasewich E.R. 1995, Wave propagation
tain narrowband signals, concentrated near resonant in heterogeneous porous media: a velocity-stress, finite
frequencies with higher velocity. difference method. Geophysics, 1995, vol. 60, p.327340.
Landau L., Livshits E. 1988, Theoretical Physics, vol. 6:
Hydrodynamics. (in Russian).
Nikolayevsky V.N. 1984, The mechanics of porous and
6 DISCUSSION fractured media M.: Nedra, 232p. (in Russian).
Rapoport M.B., Rapoport L.I., Ryjkov, V.I. 1992, Usage
Qualitative differences between the pulses reflected of seismic wave absorption method in exploration of
from the top and bottom bounds of the PL is explained hydrocarbons: Abstract of papers, 54 EAEG Meeting,
by frequency dependence of refraction indices at the Paris.

272
Copyright 2005 Taylor & Francis Group plc, London, UK
Rapoport, M.B., Rapoport L.I., Ryjkov, V.I., Parnikel V.E., Shermergor T.D. 1977, The elasticity theory of the micro-
Kately V.A. 1994, Method AVD (Absorption and Velocity inhomogeneous media. M.: Nauka 400p. (in Russian)
Dispersion): Testing and Using on the oil deposit in West- Shapiro S.A., Hurbal P. 1999, Elastic Waves in Random
ern Siberia, Abstract of papers, 56 EAEG Meeting, Vena. Media., Springer-Verlag Berlin Heidelberg. 187 p.
Ryjkov, V.I., Rapoport, M.B. 1994, Study of a seismic inelas- Tsiklauri D., Beresnev I. 2001, Enhancement in the dynamic
ticity from VSP: Abstract of papers, 56 EAEG Meeting, response of a viscoelastic fluid flowing through a lon-
Vena. gitudinally vibrating tube Phys. Rev. E 2001, vol. 64
Rapoport, M.B., and Ryjkov, V.I. 1994, Seismic velocity 036303
dispersion: An indicator of hydrocarbons: Abstracts of Yamamoto T., Turgut A. 1988, J. Acoust. Soc. Am., vol. 83,
papers, 64 SEG Meeting, Los Angeles. 1744 p.
Rapoport, M.B., Ryjkov, V.I., Rapoport, L.I., Girshgorn, Zhi-Wen Cui, Jin-Xia Liu, Ke-Xie Wang 2003. Elas-
L.Sh., etc. 1995, Interpretation of seismic inelasticity tic waves in non-Newtonian (Maxwell) fluid-saturated
effects in oil and gas prospecting., 65th Ann. Internt. porous media. Waves in Random Media 2003, v. 13,
Mtg., SEG. p. 191203

273
Copyright 2005 Taylor & Francis Group plc, London, UK