Transport in Porous Media 4 (1989), 147-184.

147
9 1989 by Kluwer Academic Publishers.

Review Article:

Transient Flow of Slightly Compressible Fluids

Through Double-Porosity, Double-Permeability
Systems- A State-of-the-Art Review
Z.-X. CHEN
Research Institute of Petroleum Exploration and Development, Beijing, People's Republic of China

(Received: 16 June 1988)

Abstract. The theory of transient flow of slightly compressible fluids through naturally fractured
reservoirs based on the double porosity conceptualization is summarized. The main achievements in
the theory of fluid flow in leaky aquifer systems which are closely related with the double-porosity,
double-permeability problems are also addressed. The main emphasis of this review is the analytical
treatment of these problems.

Key words. Double-porosity, double-permeability, Barenblatt-Zheltov, Warren-Root, naturally

fractured reservoir, aquifer system.

1. Nomenclature
A surface area
c total system isothermal compressibility
C wellbore storage constant
D ratio of matrix system permeability to fracture system permeability
h thickness of reservoir
hI thickness of fracture
hm thickness of slab matrix block
/~ ith order modified Bessel function of the first kind
ith order Bessel function of the first kind
k permeability
Ki ith order modified Bessel function of the second kind
p pressure
q flow rate of the well
q* interporosity flow rate
r radial coordinate
rm radius of spherical matrix block
148 z.-x. CHEN

s skin factor; drawdown

S storage, S = rbch
t time
T transmissivity, T = kh/l~
x Cartesian coordinate
y Cartesian coordinate
Y~ ith Bessel function of the second kind
z vertical Cartesian coordinate; variable of Laplace transform
a geometric factor controlling the interporosity flow
3' 1.781076, Euler's constant
~/ hydraulic diffusivity
)t interporosity flow parameter
/x viscosity
v eigenvalue
th porosity
to ratio of storage capacity of the fracture system to total storage capacity

Subscripts
0 initial condition
1 matrix system
2 fracture system
D dimensionless quantity
e external boundary of the reservoir
f fracture
m matrix block
ml good matrix
rn2 poor matrix
w well

2. Introduction
Naturally fractured reservoirs occur worldwide and provide a large production of
oil and gas. The essential characteristic of naturally fractured reservoirs is that
the main storage for reservoir fluids is the matrix, while the main transport
medium is the fractures. Therefore, the behavior of naturally fractured reservoirs
is radically different from that of conventional reservoirs formed by intergranular
porosity. This means that the conventional reservoir engineering techniques
based on the classical theory (Muskat, 1937a, 1949) of fluid flow through
homogeneous porous media are insufficient in this case. Thus, to meet the
requirement for efficient development of naturally fractured reservoirs a more
complex model of flow through porous media must be considered.
TRANSIENT FLOW OF SLIGHTLYCOMPRESSIBLEFLUIDS 149
3. Discrete Approach
Initially, the approach adopted for naturally fractured reservoirs was an
enumerative one, i.e., first studying the flow behavior in an individual fracture
(e.g., Lomize, 1951; Baker, 1955; Huitt, 1956), then studying the flow behavior
in relatively simple and regular fracture networks with definite sizes and
configurations (e.g., Romm, 1966; Snow, 1969; Parsons, 1966; Wilson and
Witherspoon, 1974), and finally considering the flow behaviour in a naturally
fractured reservoir itself. This approach seems quite natural and reasonable at
first glance, and has been used up to now. For example, at the Lawrence Berkeley
Laboratory one research group has been developing this approach to model fluid
flow in naturally fractured reservoirs for field studies for several years (Long et
al., 1982; Long and Billaux, 1987).
This approach may be suited to those situations where only several fractures
are of significance. And even in such relatively simple situations, it is still not
feasible to identify the specific geometric characteristics of all significant frac-
tures in the reservoir. To try to overcome this difficulty, a stochastic approach to
modeling mass transport in fractured reservoirs has been proposed by Schwartz et
al. (1983). The discrete approach is not appropriate for investigating transient
flow behavior in a naturally fractured reservoir, where a number of matrix blocks
with different sizes and irregular shapes are separated by numerous fractures
randomly distributed throughout the reservoir. In fact, even if the mathematical
and computational difficulty in connection with solving the boundary-value
problem of unsteady state flow through a naturally fractured reservoir with highly
complex fracture networks could be overcome, the sizes and configurations of the
fracture networks themselves can hardly be sufficiently defined from the rather
limited field data which is available.

4. Empirical Engineering Approach

About 30 years ago, Pollard (1959) suggested a method to evaluate acid
treatments in fractured limestone fields from pressure buildup analysis. His
method was based on the assumptions that in a fractured reservoir, fluid flows
into the well through three regions successively: from the matrix blocks into the
fractures; in the fractures; and across the skin region near the borehole face into
the well. An expression for pressure drawdown consisting of three exponential
terms, each of which corresponds to the pressure drop in one of the three regions
respectively, was then derived. This method was used to estimate the pore
volume of fractures. Subsequently, Pollard's method was extended by Pirson and
Pirson (1961) to the calculation of the matrix pore volume.
Although the graphical technique of Pollard and the formulas of Pirson and
Pirson for the interpretation of well test data may have some success, their
approach has the defect of being too approximate and, as a consequence, is
susceptible to error (Warren and Root, 1963; Kazemi, 1969).
150 Z.-X. CHEN

5. Continuum Approach
A radically different approach for naturally fractured reservoirs appeared in 1960
when the concept of flow in two overlapping continua was proposed by Baren-
blatt and Zheltov (1960), thus founding a new theory of fluid flow through
naturally fractured reservoirs - a theory of flow through a double-porosity
medium.
It is obvious that the complexity of naturally fractured networks destine the
continuum approach as the most appropriate method for this problem. A
naturally fractured reservoir can be considered as a composite of two porous
systems: the matrix block system with high storage capacity and low permeability,
and a fracture system with low storage capacity and high permeability (Figure
1). Flow occurs through the two systems separately and an interporosity flow
takes place between them. If we imagine that the interface between fractures and
blocks becomes impermeable, the system will act just like a conventional porous
medium. Here the fractures play the role of pores and the blocks play the role of
grains. The continuum approach which was used so successfully in the classical
theory of fluid flow in a conventional porous medium may be used here as well.
As we know, in a continuum approach some average characteristics of the
medium and the flow taken over a representative elementary volume (REV) are
introduced, and the basic laws governing the process are formulated in terms of
these average characteristics.
In the case of the classical theory of fluid flow through an intergranular porous
medium the main average characteristics are porosity, permeability, pressure,
flux, and the REV should be large enough in comparison with the individual pore

FRACTURE

MATRIX

VUGS MATRIX FRACTURE

(a) (b)
Fig. 1. Schematic presentation of naturally fractured reservoirs: (a) real fractured reservoir rock
(after Warren and Root, 1963) and (b) schematic presentation of fractured reservoir rock (after
Barenblatt and Zheltov, 1960).
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 151

1.0-
l.Ev I i I
>-
IScale 11 I
I I I.Ev
03
O ] I [Scale 21
n-
O
a.

0
AVERAGING VOLUME SIZE
Fig. 2. Two scales of REVs in naturally fractured reservoirs (after Shapiro, 1987).

size to include more than enough porous and grains so that it is possible to take a
statistical average over the REV. At the same time, the R E V should be small
enough compared to the flow domain so that it can be considered a mathematical
point.
In applying the continuum approach to fluid flow through a double-porosity
medium, one needs to introduce pairs of characteristic properties in each R E V
(or at each mathematical point), one for the matrix block system and another for
the fracture system, and these should be somewhat related. That is different from
all the classical cases. The double-porosity medium is thus considered as a
composite of two continua which are overlapping and mutually communicating.
Obviously, the R E V in the double-porosity medium case is a volume at the
macroscopic level with a much larger scale than the R E V for conventional
porous media. For double-porosity media an R E V should contain numerous
matrix blocks and fissures, while for a homogeneous medium only numerous
grains and pores should be contained in an REV.
Figure 2 illustrates the two scales of REVs; Scale 1 is for the conventional
porous medium inside a matrix block, and Scale 2 is for the double porosity
medium. A detailed analysis of the discrete and continuum conceptualizations of
fractured rock can be found in Shapiro (1987).

6. Barenblatt-Zheltov Model
Using the continuum approach, neglecting the inertial effect and assuming a
pseudosteady-state interporosity flow term, a system of equations describing flow
of a slightly compressible homogeneous fluid through a double-porosity medium
was obtained by Barenblatt and Zheltov (1960) as follows:
152 z.-x. CHEN

k-2div(grad pl) - q* = 051c~ op~

Ot ' (la)
/x

k2 div(grad P2) + q* Op2

= 052c2 3t '
(lb)
/x

where subscripts 1 and 2 denote the matrix block system and the fracture system,
respectively. This system of equations is similar to that describing heat transfer in
a heterogeneous medium (Rubinstein, 1948).
For interporosity flow, a pseudosteady-state regime, i.e., the process is not
explicitly time-dependent, was assumed in view of that the exchange of fluid
between the two systems is carried out mainly under a sufficiently smooth change
of pressure. A dimensional analysis results in the following expression for the
interporosity flow

q, = a (pa _ p2), (2)

/x

where a is a characteristic of naturally fractured reservoirs.

Equations (la), (lb) and (2) represent a system of coupled linear equations.
The boundary conditions at the wellbore are also coupled if the well is producing
at a constant flow rate. The double coupling makes the problem difficult to solve.

7. M o s t Simplified Case
Immediately after Barenblatt and Zheltov proposed their model, still maintaining
the most essential characteristics of the naturally fractured reservoirs, Barenblatt
et al. (1960) considered a highly simplified case where both the flow through
the matrix block system and the storage capacity in the fracture system are
negligible. Such a medium is called a 'fracture-porous' medium. In this case,
taking Equation (2) into account, Equations (la) and (lb) reduce to
Opl
o'(P2 -- Pl) = /'Z051Cl 0 7 ' (3a)

k2 div(grad p2) + o~(pl - P2) = 0 (3b)

from which a third-order equation satisfied by both pl and p2 yields
Opi k2 0 k2
[div(grad p ~ ) ] - - - div(grad p~) = O, i = 1, 2. (4)
Ot e~ Ot 1~05~cl

When a --~ % Equation (4) reduces to the classical diffusivity equation describ-
ing flow of a slightly compressible fluid through a conventional porous medium
with k2 and 051 as its permeability and porosity.
Equation (4) possesses some special properties first revealed by Barenblatt et
al. (1960) and then clarified by Barenblatt (1963). These are that sometimes an
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 153
instantaneous jump reestablishment of pressure in the fracture system may take
place as the operating condition is changed (e.g., a rate change); and any jump in
pressure and its spatial derivative in the matrix block system does not disappear
instantaneously, as in the case of conventional porous medium, but dampens
exponentially with time. The formulation of the basic boundary-value problems
of Equation (4), which is closely related to the above-mentioned special proper-
ties, was discussed in detail by Barenblatt et al. (1960) and finally clarified by
Barenblatt (1963).
Some exact line-source solutions to Equation (4) for both linear semi-infinite
and axisymmetric infinite domains were obtained by Barenblatt et al. (1960) by
using the Laplace transformation method. Some mistakes involved in these
solutions resulting from some misunderstanding of the formulation of the boun-
dary-value problem were indicated by Barenblatt (1963) himself after a short
time. The solution for the axisymmetric case is

= + qlx ~Jo(or) 1-exp do (5)

P2 P0 ~ 2 h 19 1 -~- D2/..s

which approaches with increasing ~t/(IXq~lCl) asymptotically the Theis solution

(Theis, 1935) for a homogeneous porous medium, with k2 as its permeability and
~bl as its porosity. Thus, one of the most important characteristics of the fluid flow
through a naturally fractured reservoir is that there is some lag in the transient
process compared with the case of a conventional porous medium. The charac-
teristic time of the lag is I x l c l / a . A method to estimate this important parameter
using well testing data was proposed by Ban (1961).
An exact line-source solution for a finite naturally fractured reservoir with
kl = q52 = 0 was given by Avagan (1967a), where several approximate analytical
solutions were also given and compared numerically with existing corresponding

Fig. 3. Schematic diagram of a bounded, closed reservoir with a number of wells (plan view).
154 Z.-X. CHEN

exact ones. The same problem was solved by Chen (1983) for a more general
case when a finite reservoir with any shape was penetrated by a number of wells
with finite wellbore radius (Figure 3). In the same paper (Chen, 1983), an exact
solution for the axisymmetric case was obtained, and a method of estimating all
the reservoir parameters by a number of well tests was presented.
It is worthwhile noting that in order to formulate the problem for the most
simplified case of kl = 0, r = 0, it is necessary and sufficient to impose an initial
condition on the matrix block systems and to impose boundary conditions on the
fracture systems, as was done by Avagan (1967a) and Chen (1983). Obviously, to
impose an initial condition on the fracture system is overdefining the problem
both physically and mathematically, and the instantaneous jump reestablishment
of pressure in the fracture system is a result of the improper mathematical
treatment of the problem.

8. Warren-Root Case
In 1963, an idealized model of a naturally fractured reservoir was developed by
Warren and Root (1963). Unsteady-state flow in this model was also investigated
in this paper. The model presents the naturally fractured reservoir as an idealized
system formed by identical rectangular parallelepipeds separated by an orthog-
onal network of fractures. The flow is considered to take place in the fracture
network, and the parallelepipeds feed the fracture network under a pseudosteady-
state flow condition. The system of equations resulting from this idealized model
is just a special case of Equations (la) and (lb), where kl is taken as zero

OP-A (6a)
a ( p 2 - Pt) =/~r Ot '

-- 3P2
k2 div(grad P2) + a(pl - p2) =/xq~2c2 ~-. (6b)

The significance of the idealization made by Warren and Root is that the
parameter a can be specified directly from the matrix permeability, size and
shape of the blocks. However, as Van Oolf-Racht (1982) has point out, it should
not be understood that the Warren-Root solution is limited to regularly shaped
blocks with the same size, as is sometimes mentioned in the literature, as a weak
point of the result. In fact, in order to establish such a mathematical model, the
assumption of any regular shape and identity of matrix blocks or any regular
pattern of fractures is not necessary and the only requirement is that demanded
by the continuum approach, as mentioned above. The assumed identical size and
rectangular shape of the matrix blocks are no other than some average charac-
teristics of the geometry of naturally fractured rocks.
Introducing dimensionless parameters

2 "rrk2h(po - Pi)
pD, = (i = 1, 2), (7a)
q/~
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLEFLUIDS 155
r
rD = - - , (7b)
rw

k2t
tD = tx(49,cl + 492c2)r2, (7c)

~2 C2
o) - (7d)
4,1c1 + 4~2c2'
ar~
h = (7e)
k2
Equations (6a) and (6b) for the axisymmetric case can be rewritten in dimension-
less form

a (pl)2 -- PD1) = (1 - o)) OpD1 (8a)

OtD'
1 0 (rD Opm'~ OpD2 (8b)
rD OrD O r b / + A ( p D I - p D 2 ) = w OtD
which was approximately solved with the following initial and boundary con-
ditions by Warren and Root by means of the Laplace transformation method.

pD2(~176
tD) = 0 (8C)
for the infinite reservoir or

c)pD2 [ =0 (8c')
Orb I roe

for the finite reservoir

apm [ =
(Sd)
0rD ] ,D=I
pDl(ro, O) = pD2(ro, O) = O. (8e)

For the infinite reservoir, the solution in the Laplace space is

Ko( rDg[ zf ( z )])

(9)
PD2 = z 4[ zf( z)] KI (,/[ zf( z)])'
where

~o(1 - o J ) z + A
f(z) - (1 - , o ) z + a

The inversion of Equation (9) is difficult to obtain. By considering only the first
terms of the ascending series of the Bessel functions K0 and K1, an asymptotic
analytical solution in physical space was obtained by Warren and Root (1963) as
156 z.-x. CHEN

pDz = 89{ln tD + O.80908 + Ei [ AtD AtD

where

-Ei(-x) = f~ exp(-__u) du.

~x u

Therefore, it has been found that two additional parameters, the ratio of storage
capacity of the fracture system to total storage capacity, to, and the interporosity
flow parameter, A, are sufficient to characterize the behavior of a naturally
fractured reservoir. The pressure drop or buildup curve in a semilog plot has two
parallel straight-line segments whose slopes are related to k2, and the vertical
separation of the two segments are related to to (see Figure 4). Warren and Root
(1963) presented a technique for analyzing buildup data to estimate the
parameters, k2, to, and A.
However, Odeh (1965) noted two years later that for all practical purposes,
one cannot distinguish between naturally fractured and homogeneous reservoirs
from pressure buildup or drawdown plots. It turned out that he based his

&
o
o
u

~,,~

Q.

o
z

310
DIMENSIONLESS TIM&, I D

Fig. 4. Warren-Root solution for an infinite reservoir and some particular values of to and A (after
Warren and Root, 1963).
TRANSIENT FLOW OF SLIGHTLYCOMPRESSIBLEFLUIDS 157
conclusion on calculations using parameters of some particular naturally frac-
tured reservoirs. The asymptotic solution of Warren and Root and the technique
suggested by them for analyzing the buildup data to evaluate the desired
parameters, were adopted by industry as a foundation for well-test analysis of
naturally fractured reservoirs (Matthews and Russell, 1967; Earlougher, 1977).
In 1969, an exact line-source solution of the Warren-Root case for an infinite
reservoir was obtained by Kazemi et al. (1969) by use of the Laplace trans-
formation method. In the same year, removing the assumption of pseudosteady-
state flow from the matrix to fractures, Kazemi (1969) also developed an ideal
theoretical model of a naturally fractured reservoir with a uniform fracture
distribution which we will discuss later. In these papers, Kazemi (1969) and
Kazemi et al. (1969) extended the Warren-Root solution to interpret interference
test results, and indicated how to separately estimate all key parameters, includ-
ing kl, k2, 4h, ~b2 and ~, by combining buildup test results with results of the
interference test. Some methods for estimating the interporosity flow parameter,
)t, were presented and subsequently improved by Uldrich and Ershaghi (1979),
Bourdet and Gringarten (1980), and You and Chen (1986). The first semilog
straight-line segment can exist only at very early times and is usually obscured by
wellbore storage. In the case when the first straight-line segment is not apparent,
it is difficult to use the Warren-Root analysis technique.
In the work of Crawford et al. (1976), the Warren-Root solution was combined
with a nonlinear, least-squares regression technique to analyze field buildup data.
A significant work on the analysis of pressure data influenced by wellbore
storage and skin effects, was presented by Bourdet and Gringarten (1980). This
paper contains a set of type curves for identification of flow periods and
estimation of reservoir parameters. Bourdet and Gringarten showed that double-
porosity behavior (PD) is controlled by the independent variables tD/CD, Co e 2s,
o) and )t e -2s, and that it is possible to represent the behavior of a well with well-
bore storage and skin in an infinite reservoir with double-porosity behavior as a
combination of the homogeneous behavior of each constitutive porous medium
with wellbore storage and skin at the well and the behavior during interporosity
flow from the matrix block system into the fracture system (see Figure 5). A
successful application of the new type curves to a field case was reported by
Gringarten et al. (1981).
The double-porosity behavior may often be analyzed using a homogeneous
model with appropriate boundary conditions. An efficient way to distinguish
between homogeneous and heterogeneous (including double-porosity) behavior
is to examine a log-log plot of the derivative of po (Figure 6). Such a new set of
type curves for interpreting well tests in naturally fractured reservoirs was
introduced by Bourdet et al. (1983).
In 1979, Mavor and Cinco-Ley (1979) presented their comprehensive study on
double-porosity systems, taking wellbore storage and Skin effects into con-
sideration. After the corresponding solutions in Laplace space were obtained,
158 z . - x . CHEN

102 ' -' ' ' '

C^e2S - . . . . . 103o
~ =_e_____ .,~"~'-~"l . . . . . . . . ~. . . . . . . . 10"'~
APPROXIMATE ~ J ' ," ---IOIs
. . . . SiAm OF SE -LOG / ~ - J --- ...... ~o

<,.10 _,o,
.,W t'~
// ~ ~ :ii ~
43 L " " 10-4
. ' . . . . . . . . . . . . . . 10"2

. . . . . . . . . . . . .

10"1 1 10 102 103 104

*._o.o o.oooles ~
CD #J C

Fig. 5. Wellbore storage and skin type-curve in a double porosity reservoir (pseudosteady-state
interporosity flow) (after Bourdet and Gringarten, 1980).

they were numerically inverted into real space by means of a numerical Laplace
invertor, proposed by Stehfest (1970). Subsequently, Da P r a t e t at. (1981)
applied the solution for constant bottomhole pressure given by Mavor and
Cinco-Ley (1979) to carry out decline curve analysis using type-curves for
double-porosity reservoirs.

lO i i

HOMOGENEOUS BEHAVIOR
"h I
O C0=i s . iOe "~
I=l 0.5
APPROXIMATE START /
OF I N F I N I T E ACTING
HOMOGENEOUS BEHAVIOR
O DOUBLE POROSITY
a.
'i= I0"I BEHAVIOR (Coe 2$}~ i 0 I
(CDr t $ l f . ~ , , i . 5 1 0 5
X o - z s . 510 -ao

i i i i
IO'|
I0 I0 i 103 104 108
tDIC o
Fig. 6. Derivatives for homogeneous and double-porosity behavior (after Oringarten, 1984).
TRANSIENT FLOW OF SLIGHTLYCOMPRESSIBLEFLUIDS 159
In 1977, overcoming the difficulties encountered in the inversion of the
Laplace transformation, Jiang (1977) obtained the exact solutions of the Warren
and Root cases to all the basic boundary-value problems for both infinite and
finite axisymmetric domains with a finite internal radius. Jiang successfully
arrived at the calculation of the following integrals

1 f rD+i~&tKo(ro~x)
dz
2 rri JrD_i~ xTxxK~(~x)
and

1 ~ rD+im eZta1/o.l(rD,roe,~X)
Zt dz '
27ri % i~ ~Xatrl,l(1, rDe, ~XX)

where

~,.,. (t~,/3, y) = K.,(c~y) In (/3y) + (-1)m+"+llm(c~y)K.(/3y)

and
z(1 +az)
X--
1 +bz
As we see, the integrands are so complicated that analyzing and calculating the
number and distribution of their singular points on the complex plane become
very difficult. Usually, the complex integrals can be calculated by choosing a
contour such that only a finite number of singular points are inside the contour
and Cauchy's integral formula can be used. In Jiang's case, however, there will
always be an infinite number of singular point inside any allowable contour and,
therefore, a new method must be developed to solve the problem.
The solutions obtained by Jiang are

pm(rD, to) =--2 H(tD, y)F(rD, y) dy (11)

7r
for an infinite reservoir;

pD2(rD, to) = 2 [ tD + (b - a)[1 -exp(-tD/a)]

r2De- 1

+ ~ H(to,ak)Gl(ro, rDe,Otk)] (12)

k=l

for a finite reservoir with no-flow outer boundary; and

pm(rD, to) = 2 ~ H(tD, flk)G2(rD, rDe, ilk) (13)

k=l
160 z.-x. CHEN
for a finite reservoir with a constant-pressure outer boundary; where

1--0)
b- , a = rob,
h

H(to y ) = ~, A'(Y) II-exp[-o'~(y)to]},

i=1 ~ / ~ ~-

ni(y)=~a
1 [b-(-1)' b(l+by2)-2a
x/[(~Tl_-~y-~-~aay2] j
] ( / = 1,2),

o's(y) = ~ a { 1 + b y 2 - (-1)~/[(1 + by2)2-4ay2]} (i = 1, 2),

G~(rD, rDe, c~k)= @o,~(rD, rD,, '~k)

rDel~JO,l(rDe, 1, ak)-- (I)0,1(1 , rDe , Olk)'

G2(ro, ro., elk) - @o,o(ro, ro., fig)

roe@~A(1, rD., i l k ) - @o,o(1, ro., ilk)'
@m,,(a, /3, y) = Ym(ay)J,([3y) - Jm(ay) Y. (/3y),

= @0,1(rD, 1, y)
F(rD, y) j ~ ( y ) + y~(y)
2

ak is the kth positive root of the equation

@1,1(1, rD~, y) = 0

and ~k is the kth positive root of the equation

@o,1(1, rDe, y) = 0.
The dimensionless bottomhole pressure subsequently calculated by Zhu et al.
(1981) according to Jiang's exact solution, shows that the Warren-Root asymp-
totic solution is in excellent agreement with the exact one, except at very early
times (see Figure 7).
The exact solution of the Warren-Root case for finite reservoirs with a no-flow
outer boundary was also obtained later by Chen and Jiang (1980) using the
method of the separation of variables in a muGh easier way. In this work, the
Warren-Root model was also used to study a mor,e general case of flow in a finite
reservoir with a no-flow outer boundary of any shape and with numerous wells
arbitrarily located (Figure 3). Based on the structure of the solution obtained by
Chert and Jiang (1980) for the general case, the average pressure of both the
fracture and matrix block systems, and the rate of crossflow between the two
media for such a general geometry of reservoir, were obtained by Chen (1982).
The exact solution of the Warren-Root case, including wellbore storage and skin
effect, was obtained by Luan (1981a) in terms of the Laplace transformation
TRANSIENT FLOW OF SLIGHTLYCOMPRESSIBLEFLUIDS 161

PD

10

Exact solution

_ _ _ Approximate solution

A= 10 - 5

#" f -
-4 / \v~ 2 0 2 4 6 8 Igto

Fig. 7. Comparisonof Warren-Root approximate solution with Jiang's exact solution for an infinite
reservoir (after Zhu et al., 1981).

following the method in Jiang's solution (1977), and later by Liu et al. (1987)
using the method of the separation of variables.
Under the stimulus of Jiang's work (1977), Several other analytical solutions of
fluid flow in double-porosity media, including a solution of flow towards a
partially penetrating well by Zhang (1982), a solution of flow through a com-
posite reservoir by Wu and Ge (1981) and a solution for non-Newtonian fluid
flow by Luan (1981b), were obtained in China in the early 1980s. The problem of
flow to a partially penetrating well including weUbore storage, skin effect and a
nonequilibrium initial condition in the wellbore, was solved by Dougherty and
Babu (1984) using the Stehfest algorithm.
Some attention has been paid to the problem of identifying the parameters in a
naturally fractured reservoir through field data. In 1978, a theoretical mode] for
flow in a bounded, closed naturally fractured reservoir was proposed and dis-
cussed by Sun and Bai (1978) using spectroscopic eigenvalue analysis. The
optimization method to automatically identify reservoir parameters developed for
conventional reservoirs, was generalized to double-porosity systems by Lang and
Li (1981). A technique of parameter estimation using cybernetics was presented
by Qi and Zhang (1986).
162 Z.-X. CHEN
9. Double-Porosity, Double-Permeability Case
The double-porosity model which neglects flow within the matrix block system,
generally yields satisfactory results because the matrix permeability is usually
much less than that of the fracture system in a naturally fractured reservoir.
However, in order to estimate the limits of validity of solutions based on the
double-porosity model and to study the behavior of a naturally fractured reser-
voir when the contrast between the two permeabilities is not significant, it is
necessary to solve the original model proposed by Barenblatt and Zheltov
[Equations (la) and (lb)].
The Barenblatt-Zheltov model in its complete form was first studied by
Bondarev and Nicolaevsky (1966). The conditions under which Equations (la)
and (lb) can be reduced to Equations (3a) and (3b) were investigated, and an
approximate analytical solution of Equations (la) and (lb) for a semi-infinite
linear reservoir with constant pressure at the inlet was presented. An approximate
analytical solution of Equations (la) and (lb) for an infinite reservoir with a well
producing at a constant rate was obtained by Avagan (1967b). Some known
methods of well test analysis were generalized to estimate the parameters of
double-porosity, double-permeability media (Shalimov, 1966, and Kutliarov,
1967).
In 1980, the Barenblatt-Zheltov model was first solved rigorously by Chen and
Jiang (1980) by proper decomposition of the problem and using the method of
the separation of variables. A bounded, closed naturally fractured reservoir of an
arbitrary shape with spatially varying properties and a number of wells arbitrarily
located and produced or injected at given flow rates is considered (Figure 3).
The external boundary of the reservoir is denoted by OfL , the boundary of
the ith well is denoted by 0ft~ ) and the domain between 0fL and 0 ~ )
(i = 1, 2 , . . . , n) is denoted by ~. The problem is then formulated as follows

~x(klf~Pl~-[-L(klOPl~-Jrol(p2-Pl)~-/.s~ t 1,
Ox] Oy Oy/
(x, y) E ~~, t~0,
(14a)

0 k2 k2 +a(pl-pz)=tX4)zCz (x,y)et2, t>0,

Ox ~x ] Oy Oy ]
(14b)

OPl I = OP--~2I =0, (14c)

art of~+ cqn 0~+

( i = 1 , 2 , . . . , n), (14d)
m" On OnJ - ,~ q(i)(t)'

pl = p2 = x(i)(t), (x, y) e 0~(_/),

X(0(t) is to be determined, (i = 1, 2 . . . . . n), (14e)

TRANSIENT FLOWOF SLIGHTLYCOMPRESSIBLEFLUIDS 163

pl(x, y, o) = plo(x, y), (14f)

p2(x, y, O) = p2o(X, y), (14g)

where q(i)(t) is the flow rate of the ith well.
This is the Barenblatt-Zheltov's system of equations with equipotential surface
boundary conditions. As we see, both the system of equations and the boundary
conditions are coupled.
A structure of the solution for the general case was obtained, resulting in a new
technique to estimate oil reserves in a naturally fractured reservoir. For a case of
a single well located at the center of a circular drainage area, an exact solution
was obtained by constructing a set of special quasi-Bessel functions. The
numerical calculation of this solution is not given in this paper, because the
quasi-Bessel functions were not well understood at that time.
Liu and An (1982) solved the Barenblatt-Zheltov model numerically. Un-
fortunately, their results did not reveal the important early-time pressure res-
ponse. The Barenblatt-Zheltov model, including wellbore storage and skin effect,
was investigated by Bourdet (1985) resorting to the Stehfest algorithm of
numerical inversion of Laplace transforms. In 1987, Chen and You (1987)
adopted a new approach to avoid the quasi-Bessel functions. The approach was
based on first solving the system of equations for an arbitrarily varying pressure,
pw(t), as its inner boundary condition, and subsequently pw(t) was determined to
satisfy the constant flow rate condition as given in the original problem, resulting
in an integral equation containing only ordinary Bessel functions. This integral
equation can be solved numerically, thus allowing the studying of flow behavior
in naturally fractured reservoirs including fluid flow through the matrix block
system. The exact solution of double-porosity, double-permeability media
expressed in terms of ordinary Bessel functions, was first indicated by Lin and
Chen (1987) as a special case in a still unpublished paper where a multiple porous
medium problem was dealt with and an orthogonal matrix transformation was
proposed to make the coupled system of n equations uncoupled. This is further
discussed below.
In 1987, some exact solutions of the Barenblatt-Zheltov model including
wellbore storage and skin effect were obtained by Liu et al. (1987) for a centered
well producing at a constant rate in a reservoir with no-flow or constant-pressure
outer boundary, using the method of separation of variables. When the wellbore
storage was included, the resulting eigenvalue problem became very unusual.
One of the boundary conditions contains eigenvalue itself. Some new operators
were introduced to eliminate this problem. Similar to Liu and Chen (1987), a
Jordan transformation was used to make the coupled system of equations into a
system of equations consisting of two independent Bessel equations of zero order
and, as a consequence, the solutions are expressed in terms of ordinary Bessel
functions. For the no-flow outer boundary case, the solution is expressed as
follows
164 Z.-X. CHEN

P m = O[ tD -t r2D4 r2~ln rD +I-~-[G~ In rD~-

3r4e-- 2r e-- '] 0

16 - ~--~ (I - ~2) [ O2(1) + ~[ --CDOf~--

-- E Hlj[~ l( rD,/~11) q- 1~2jBIIOPO,I(rD,

* P~2j)] e-~'D, (15a)
]=1
_tr~ r~eln ] D [r,~D~
PD2 = O[ tD 4 rD - i ~ - ~ [ G o - GI(rD)]+ 02 In r D , -

3r4, - 2r2, - 1] 0 [ + Go] _ cDof _

16 - ]-~ (1 - r G2(1) 0~1

-- E Hli[txljAud~O,l(rD,
~ /&li) "at-~2B,iOO, I(rD, /x2i)] e --u.t,D, (15b)
j=l

where

2~r(kl + k2)h
p~, - (po-p~) ( i = 1,2),
qt*
kl q- k2
tD -- /~(r + 4)2c2)r 2 l,
,~r2 k~ C
~.- D- C~=
kl + k2' k2' 2 ~r( (/~1C1-1-(~2c2) r2.

[Other symbols are defined in Liu et al. (1987).]

It has been shown that the quasi-Bessel functions constructed by Chen and
Jiang (1980) are linear combinations of ordinary Bessel functions, and the exact
solution obtained in Chen and Jiang (1980) is identical to the solution obtained
here for the reservoir with no-flow outer boundary for the special case of
neglecting wellbore storage and skin effect.
Numerical calculation of these exact solutions were made and the behavior of
the double-porosity, double-permeability systems was studied as a function of
various reservoir parameters. The limitations of the double-porosity model was
also examined. The permeability ratio D weakens the effect of double-porosity
media revealed by the Warren-Root model. The early straight-line segment of
semilog dimensionless bottomhole pressure curves moves closer to the late
straight-line segment and, in general, they do not become parallel to each other
as in the Warren-Root case (see Figure 8). The conditions for using the
well-known method to estimate oJ based on the Warren-Root model, were
established in Liu et al. (1987). It can be used if two conditions are satisfied
simultaneously: D < 0.01 and ( 1 - w)DMo~ < 10 -6, but even if these conditions
are satisfied, the error may be up to several percent.
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 165

P~
IOF

rDe:5000
~.= I0-5
co= 0.01
S=CD= 0

i rD = 0
rD = 0.00

r%
r

-1 0 1 2 3 4 5 6 7 logt~
Fig. 8. Influenceof permeabilityratio D on dimensionlessbottomholepressure drop (afterLiu et al.,
1987).

All the existing exact solutions, such as the solutions of Hurst (1934) and
Muskat (1934) for homogeneous porous media, that of Jiang (1977) for double-
porosity media, and the Chen-Jiang solution (1980) for double-porosity, double-
permeability media, are special cases of the general solutions.
The exact solution of the Barenblatt-Zheltov model for a constant pressure
condition at the bottomhole, was obtained by Ge and Wu (1982) and Chen and
You (1987). This kind of problem is considerably easier to solve than the
problem with a constant flow rate at the wellbore, because the boundary
conditions are not coupled. When the skin of each porous system is included into
the inner boundary conditions, the problem becomes much more complicated.
Some exact solutions were obtained recently by Chen (1988) using practically the
same method as in Liu et al. (1987).
166 z.-x. CHEN
10. Transient Interporosity Flow Model
All the works mentioned above were based on the assumption of pseudosteady-
state interporosity flow. That means that the intensity of interflow from the matrix
block system to the fracture system is in direct proportion to the pressure
difference between the two porous systems. Obviously, this is an approximation.
The flow from the matrix block into the fracture must be a transient phenomenon
and only in due course reaches a pseudosteady-state condition. Thus, instead of
Equation (2), the interporosity flow rate should be expressed as
k,,
q* - A,,(grad Apm)int, (15)
/x
where subscript int denoting the gradient of Ap,, is taken at the interface between
the fracture and matrix block systems.
The first paper considering transient interporosity flow was given by Kazemi
(1969). An ideal model consisting of a set of horizontal fractures and a set of
uniformly spaced horizontal matrix layers with the set of fractures as the spacers,
was developed to represent a naturally fractured reservoir. The solution was
obtained numerically using an iterative alternating direction implicit procedure.
Kazemi concluded that the major results of Warren and Root are acceptable. The
model Kazemi studied numerically was general, considering both the transient
interporosity flow and flow through the matrix block system. Unfortunately,
Kazemi did not pay enough attention to the considerable divergence of the
transitional segment caused by the different regimes of interporosity flow. Fur-
thermore, he did not try to determine the limitation of the double-porosity model
caused by neglecting fluid flow within the matrix block system.
Dontsov and Boyrchuk (1971) were the first to solve the transient interporosity
flow model analytically, but, like other Soviet authors, only for the most sim-
plified case when both the flow through the matrix block system and the storage
capacity in the fracture system are neglected. The Laplace transformation was
used and the inversion was carried out approximately. Their solution showed the
existence of a half-slope straight line preceding the late straight line when plotted
on a semilog paper.
An impetus to the use of transient interporosity flow was given mainly by de
Swaan (1976). The storage capacity of the fracture system was taken into
consideration but the flow through the matrix block system was still neglected.
Two kinds of identical uniform matrix blocks, i.e., slab and sphere, were assumed.
The interporosity flow rate is expressed as
2 ft OAp2
q*(p2, t) = - a , ~ Jo - ~ - q . ( t - ~-) d~-, (16)

where
k,.
q. = - - - A,.(grad Ap.,.)i.t, (17)
T R A N S I E N T F L O W OF S L I O H T L Y C O M P R E S S I B L E F L U I D S 167

(27r+ l)7rz
Apu,,,(z, t)= 1 - ~4~=o (-1)" e_n..(z.+l)2.C~/h~cos
= 2n+-----1 h,. (18a)

for slabs; and

+ 2r,, ~__ ( - )1 n -.7m.2.,2t/. ms i n nvrr
Ap,,m(r, t) = 1 --e -- (18b)
7rr =1 n rm

for spheres.
An approximate large-time solution was obtained in the same form as the
well-known approximate solution of the radial infinite homogeneous reservoir,
but with a modified hydraulic dittusivity constant:

_ q/x In [ 4~ ]

where
1
- 1 k,,,h,,, (20a)
-4
~r klhfq-~
for the slab model; and
1
n- 1 2kmrm (20b)
-4
~qf 3ktlq'qm
for the spherical model.
No analytical description of the transition segment was given, while the
approximate solution of a radial infinite homogeneous reservoir with the same
parameters as in the fracture system was taken as the early-time solution.
De Swaan's theory exclusively involves flow properties and dimensions of the
fractures and the matrix blocks; and no extra adjusting parameters are needed for
predicting the behavior of a reservoir with known properties.
Shortly after de Swaan's paper, Duguid and Lee (1977) derived equations
governing the flow of fluid through naturally fractured media including an
acceleration term for flow through the fracture. A transient interporosity flow
term was obtained assuming slab-shaped matrix blocks. A numerical solution was
achieved using the Galerkin finite element method for a confined naturally
fractured leaky acquifer. In the same year, Boulton and Streltsova (1977b)
proposed an identical model to the slab case in de Swaan's model. The new
contribution of Boulton and Streltsova was to present an exact line-source
solution for a radial infinite reservoir as follows
c~
_ q r
sm 2,a.hrk, fo X J o [ ~ X ] [ ~ gr,~*,~]dx, (21a)

q oo r
(21b)
168 z.-x. CHEN

where
1 - exp(-/3~rtqt/h~ )
~ , . - [ ~ " h " / ( ~rhf )]/3~ + 0.5[ h,.k,./hfkt]/3,.(tan /3,. +/3~ see 2/3,.)'
9 " = cos(/3,~z/h,,.) + tan/3,, sin(/3mz/h"),

/3" is a positive root of equation

]
tan/3,.] = (hma) 2
W
(a is the parameter of the Hankel transform).
McNabb (1978) presented a fruitcake model which is similar to de Swaan's
model and is differentiated from the latter only by a more general shape of the
matrix blocks. McNabb gave an early-time point-source solution for an infinite
spherical reservoir. About the same time, some approximate analytical solutions
of de Swaan's model, expressed in the form of Theis' line-source solution with
some modified diffusivities which were quite convenient for practical use, were
given by Najurieta (1980). The features of the transition segment under a
transient interporosity flow regime were indicated. The well testing method
proposed by Najurieta makes it possible to estimate all four parameters, T~, Sin,
St and ~- [~-= h~/(4y~m)], which fully describe the reservoir behavior according
to the solutions.
In 1983, a number of authors (Streltsova, 1983; Serra et al., 1983; Cinco-Ley
and Samaniego, 1982) made an important finding that, just as Dontsov and
Boyrchuk (1971) did in their case (D = w = 0), that, instead of the plateau
transition segment indicated by Warren-Root's solution, there was a linear
transition segment whose slope is equal to half the slope of the classical parallel
semilog straight-line segments (Figure 9). These authors derived some ap-
proximate analytical solutions based on the transient interporosity flow regime.
They accordingly presented some new, useful transient well-testing analysis
methods. In most tests carried out in naturally fractured reservoirs, the first
semilog straight-line segment can hardly be obtained. It terminates at very sho(t
times and is usually obscured by the wellbore storage and skin effect. This makes
it difficult to estimate the reservoir parameters using conventional semilog
analysis methods based on the Warren-Root solution. Now, however, there are
new methods which make this much easier. When two or more of these three
semilog straight-line segments are available, a virtually complete analysis of
pressure data is possible. It should be mentioned that the existence of the semilog
straight-line segment during the transition period was first pointed out by
Bourdet and Gringarten (1980). A little bit later, a similar model with matrix
blocks of identical cubic geometry was presented by Lai et al. (1983). Solutions
were obtained, including wellbore storage and skin effects in the Laplace domain,
and were inverted numerically. The half-slope segment was also observed for
values of r smaller than 0.1 and followed by a brief segment with a slope of 2/3.
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 169

12

S,~/SI ; 100

0
101 10' 1.0~ 10 3 101 101 10 i 10'
4t D = 4rjt/r ~
Fig. 9. Comparison of pressure patterns from the pseudosteady-state interporosity flow and the
transient interporosity flow models (after Streltsova, 1983).

Moench (1984) introduced a concept of fracture skin, namely, a thin skin of

low-permeability materials deposited on the surfaces of the matrix blocks, that
impedes the interporosity flow. Solving this problem by Stehfest's numerical
Laplace invertor, Moench arrived at the conclusion that fracture skin provides a
theoretical justification for the pseudosteady-state interporosity flow ap-
proximation. Transient interporosity flow models were extended by McGuinness
(1986) to matrix blocks having sizes varying according to exponential dis-
tributions. The existence of a half-slope segment was shown using an analytical
solution in Laplace space and conditions for observing this were derived.
However, there is still some divergence of views about the semilog half-slope
segment. Ershaghi and Aflaki (1985) showed, based on both theoretical analysis
and published data, that a 1:2 slope ratio was not necessarily indicative of the
transient interporosity flow regime. It seems that, in general, the Warren-Root
model reflects the essential features of naturally fractured reservoirs more often
than not, but many aspects, including the interporosity flow regime, remain to be
studied in more detail.
Some approximate solutions for a bounded cylindrical reservoir with a central
well produced at a constant rate, and that for either an infinite or a bounded
reservoir produced at a constant pressure, in addition to welt-test analysis
methods based on them, were presented by Chen et al. (1985) and Ozkan et al.
(1987), respectively.
Type curves for the analysis of drawdown and buildup data, as well as pressure
transient data from observation wells using the transient interporosity flow model
and their applications to field data, can be found in Bourdet and Gringarten
(1980), Gringarten (1984), Deruyck et al. (1982), Lai et al. (1983), Chen et al.
170 z.-x. CHEN

(1984), and Reynolds et al. (1985). These curves are valuable but must be used
with care to avoid potential errors caused by multiple matches. Some ap-
proximate engineering methods were developed by Aguilera (1987a, 1987b) to
analyze well test data under more realistic reservoir and operation conditions.
Almost all authors assume identical, uniform and regular shapes of matrix
blocks separated by some simple regular fracture networks to study the reservoir
behavior under a transient interporosity flow regime. Such an idealization is very
helpful as a means of visualizing the problem, but it must be kept in mind that the
results deduced from such an idealization are not only applicable to the parti-
cular model used. From the point of view of the mechanics of continuum, the
idealized geometry of matrix blocks and fractures are also just average charac-
teristics over a representative elementary volume containing many matrix blocks
and fractures of all kinds. This fact has unfortunately been overlooked by some
reservoir engineers, and sometimes even by authors.
Except for the Boulton-Streltsova's solution (Boulton and Streltsova, 1977b),
which is a line-source exact solution, all the existing solutions for the transient
interporosity flow model were obtained in Laplace space and then numerically
inverted or approximately inverted into real space. Of course, such solutions
cannot be considered as being exact and, moreover, no solution has been
obtained for the double-porosity, double-permeability case assuming a transient
interporosity flow regime.

11. Triple-Porosity Model

Some attention has been paid to triple-porosity or, more general, multiple-
porosity systems. This was prompted by an investigation of the rock type often
encountered in carbonate reservoirs and observations of actual well tests in such
reservoirs that show anomalous behavior that could not be explained by double-
porosity models.
A triple-porosity model describing the flow behavior in a fissured aquifer,
where two basic rocks having good and poor petrophysical properties respec-
tively coexist, was generalized from the Warren-Root model by Closmann
(1975). The mathematical model is
Op,~l
a,,,l ( p f - p,,,1) = aq5,,,l c,,~l -~ , (22a)

0p,,,22 (22b)
amz(pf-- pm2)= OLI])m,2Crnz Ot '

- Opf
r02ps + 10ps ] + a,,,, (p,,,, - Pr) + ~,,2(p,,2 - Ps) = #xOscs ~-~" (22c)
k'L ; 0,J

Exact solutions obtained by a finite Hankel or Fourier transform with respect to

distance and a Laplace transform with respect to time, were presented for both
bounded, closed radial and linear reservoirs with constant-pressure inner boun-
TRANSIENT FLOWOF SLIGHTLYCOMPRESSIBLEFLUIDS 171

~, =0.I, toz=O.89, ~., = J64 /

I

PD

2~ A2 = ~d6

,t I I, I 9 Jt
10 7 102 103 I0~' 10s 106 107 t0 B
ID
Fig. 10. Schematics of PD vs. to for the triple-porosity model (after Wu and Ge, 1983).

dary conditions, and applied to the estimation of the water influx in a material
balance calculation problem. The exact solution for a bounded, closed radial
reservoir with constant-rate inner boundary condition was given in Liu (1981)
using the separation of variables method, while the exact solution for an infinite
case was given in Wu and Ge (1983) using the Laplace transformation method.
Instead of two parallel semilog straight-line segments for the double-porosity
systems, there are three parallel semilog straight-line segments for the case of
triple-porosity systems (Figure 10).
In 1986, a triple-porosity model with a transient interporosity flow regime was
introduced and studied by Abdassah and Ershaghi (1986). The authors included
wellbore storage and skin effect in their model and used the Stehfest's numerical
inversion algorithm to solve the problem. The application of obtained results to
interpret well-test data was discussed.
A multiple-porosity, multiple-permeability model generalized from the Baren-
blatt-Zheltov model was formulated and rigorously solved by Liu and Chen
(1987). This will be discussed in more detail below.

12. Layered Reservoirs

The behavior of a naturally fractured reservoir is analogous to the behavior of a
two-layered reservoir. A similar analogy exists between a multiple-porosity
system and a multiple-layered one. In fact, the mathematical models describing
the flow of fluids through double-porosity or multiple-porosity media are very
similar to those for two-layered or multiple-layered systems, and sometimes they
are identical. The difference is mainly in the formulation of the boundary con-
ditions. The two kinds of problems are closely related, and progress on either
problem has made it easier to solve the other. At the same time, the similarity
between the two kinds of problems makes it difficult to differentiate the naturally
fractured reservoirs from the two-layered reservoirs by well-test data alone. The
172 z.-x. CHEN
two parallel straight-line segments in a semilog plot are not only characteristic of
a double-porosity medium but characteristic of a two-layered system as well.
The theory of fluid flow through layered systems was developed much earlier
than the theory for double-porosity media, mainly in the field of hydrology. Here
the review will be limited to work that is closely related to the theory used for
double-porosity, double-permeability media.

13. Commingled Systems

The earliest mathematical analysis of oil production from a multiple-layered
reservoir appears to be the Tempelaar-Lietz work (Tempelaar-Lietz, 1961)
originally published in a 1953 Shell Oil Co. report. A bounded system consisting
of two layers, unconnected except at the well, was considered. This model is
equivalent to the Barenblatt-Zheltov model with a zero interporosity flow
parameter. Simple formulas for performance analysis were obtained based on
steady-state flow conditions.
A similar problem for n layers in an infinite reservoir was rigorously solved by
Homer in an open unpublished work using the Laplace transformation method.
The same problem was solved by Lefkovits et al. (1961) for a bounded reservoir,
and shortly afterwards by Duvaut (1961), by the same method as used by Horner.
The solution, including skin for each layer, was also presented. It has been found
that the duration of the transient period is often in orders of magnitude longer
than the transient period in a single-layered reservoir. Also, during the transient
period, the more permeable layer is depleted faster than the less permeable.
However, the reservoir will approach a pseudosteady-state and, under this flow
regime, each layer produces at a fractional rate equal to the fraction of storage
capacity it contains. Theoretical build-up curves from such commingled reser-
voirs were also examined in Lefkovits et al. (1961), showing flattening after an
initial semilog straight-line segment, then steepening, and finally flattening
towards static pressure.
Lefkovits (1961) indicates that from the buildup curve, some conventional
methods for single-layered reservoirs can also be used to estimate the average
properties of commingled reservoirs such as kh, s and static pressure, but it is not
possible to determine the properties of the individual layers. Well-test analysis for
wells producing commingled reservoirs was discussed and developed later in a
number of papers. Attention was centered on the determination of the proper
application of conventional analysis methods to the commingled-system case
(Kazemi, 1970; Cobb et al., 1972; R a g h a v a n et al., 1974; Earlougher et al.,
1974). The conclusions reached were that under well-defined conditions, the
Muskat (1937b) trial-and-error, Miller-Dyes-Hutchinson (1950) and Horner
(1951) methods can all be used to determine average parameters of a layered
reservoir, but there is no generally valid criterion for recognizing multiple-layer
systems from transient tests. This is because commingled systems do not always
TRANSIENT FLOW OF SLIGHTLYCOMPRESSIBLEFLUIDS 173
give buildup curves which are similar in shape to those mentioned above, while
some reservoirs which are not commingled may do so (Raghavan et al., 1974).
Raghavan et al. proposed a method to estimate the individual-layer per-
meabilities from the second build-up part on the Horner plot.
In 1978, a numerical Laplace invertor (Stehfest, 1970), which made many
recent advances in well-test analysis possible, was introduced into the reservoir
flow investigations by Tariq and Ramey (1978). They solved a flow problem of
multiple-layer systems of different radii without crossflow including wellbore
storage and skin effect by use of the Stehfest algorithm. It has been found that
false wellbore storage effects may appear in cases involving high permeability
contrast and a small, highly permeable layer, and that layered system data can be
analyzed under certain circumstances to yield information about the permeability
ratio and the radii of the individual layers.
Both exact and asymptotic approximate solutions for an infinite two-layered
acquifer with unequal initial pressures were presented by Papadopulos (1966).
The same problem which was extended to n layers with a no-flow outer
boundary, was rigorously solved by Larsen (1981). Methods for the estimation of
reservoir parameters from pressure transients prior to production were presented.
Based on the analytical approximation given in Larsen (1981) for the wellbore
pressure in the infinite-acting period, the same author proposed methods to
estimate both the skin and the flow capacity of each layer (Larsen, 1982).
A multilayer testing technique consisting of a number of sequential flow tests
with a production logging tool measuring the wellbore pressure and flow rate at
the top of each different layer to obtain the permeability and skin factor of each
layer for commingled layered reservoirs, was developed by Kucuk et al. (1986).

14. Systems With Crossflow

As early as 1930, de Glee (1930) theoretically analyzed the problem of the
steady-state leakage of water through less permeable layers into an aquifer that is
being pumped.
Jacob (1946) treated this problem assuming that the leakage is varying linearly
with the instantaneous difference in heads across the aquifer. This work was later
improved by Hantush and Jacob (1955) still neglecting both the storage of the
aquitard and the drawdown in the unpumped aquifer. The solution which was
obtained by Hantush (1960) included the storage of the aquitard but neglected
the drawdown in the unpumped aquifer. This was followed by a solution given by
Charny (1961) which included the drawdown in the unpumped aquifer but
neglected the storage of the aquitard. The mathematical model solved by Charny
is different from the Barenblatt-Zheltov model only in the internal boundary
conditions. For the double-porosity, double-permeability media, the internal
boundary conditions for the two porous systems must be coupled when the rate of
174 z.-x. CHEN
withdrawal is specified, while it is not necessary for two-layered media, as in the
Charny solution.
The problem of fluid flow through an infinite two-layered reservoir with
crossflow produced by a well at a constant bottomhole pressure, was solved by
Maksimov (1960). The problem was formulated separately for each homo-
geneous layer as a two-dimensional problem in the r, z space. The problems for
two layers were then linked by the interface conditions that both pressure and
fluid flux must be equal. A solution in the Laplace space was obtained by use of
the Laplace transformation in combination with a Weber transformation. A
technique to estimate the average transmissibility of the layered reservoir was
given based on the method of determining the reservoir parameters from
well-test data directly in the Laplace space proposed by Barenblatt et al. (1957).
Jacquard (1960) presented an exact solution to the problem formulated as in
Maksimov (1960) but for a well producing at a constant rate from a bounded
two-layered system with crossflow using the Laplace transformation method.
Jacquard's solution was calculated by Pottier (1961).
In 1962, three papers (Russell and Prats, 1962a; Katz and Tek, 1962; Pender-
grass and Berry, 1962) appeared dealing with bounded multiple-layered reser-
voirs with crossflow. A central well was produced at a constant pressure. In the
first two papers, the problem was formulated as in Maksimov (1960). An exact
solution for the radial case was obtained in Russell and Prats (1962a) using a
finite Hankel transformation followed by a Laplace transformation, while exact
solutions for both linear and radial cases were obtained in Katz and Tek (1962)
using the separation of variables method. In Pendergrass and Berry (1962), the
problem was formulated for the whole flow region as a two-dimensional problem.
The method of solution involved a combination of both analytical and numerical
techniques. Russell and Prats (1962b) summarized the practical aspects of these
theoretical works indicating that, except for at a very early time, a well in a
layered reservoir with crossflow behaves as in a homogeneous single-layer
reservoir with the same dimensions and pore volume as the crossflow system and
a permeability-thickness product kh equal to the total kh of the crossflow system.
Well-test data from this period can therefore be interpreted by use of homo-
geneous-reservoir theory.
In 1969, a line-source solution, including both the storage of the aquitard and
the drawdown in the unpumped aquifer, was developed by Neuman and Wither-
spoon (1969) using the Laplace transformation for an infinite radial system. The
solution obtained by them was an exact one. However, the horizontal flow in the
aquitard was neglected as usual.
In 1977, Boulton and Streltsova (1977a) studied the unsteady flow of a pumped
well in a two-layered water-bearing formation. The model considered is very
similar to the Barenblatt-Zheltov model, but the inherent pseudosteady-state
interporosity flow assumption in the Barenblatt-Zheltov model is now replaced
by a transient interporosity flow regime. However, the assumption that for each
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 175
layer the discharge per unit length is proportional to its permeability was imposed
by the authors, avoiding the real problem involved in solving such kinds of
problems. This assumption cannot be used in the double-porosity, double-
permeability case, although the line-source solution for an infinite domain
obtained by them can still be regarded as an important contribution to the
two-layered problem.
In 1987, the same problem as treated in Neuman and Witherspoon (1969) was
rigorously solved by Chen et al. (1987) without any restriction on the flow in the
aquitard. The method used was the separation of variables. A general solution for
a bounded domain of an arbitrary shape with a number of arbitrarily located wells
(see Figure 3) was presented, and an exact solution for a bounded axisymmetrical
case was obtained. It has been indicated that the usually adopted assumption of
neglecting radial flow within the aquitard introduces a significant error in the
drswdown curves only in the aquitard itself and makes the transient segment of
these curves too convex (see Figure 11).
Flow to a well in a multiaquifer system was investigated in Hunt (1985),
however, storage in the aquitards was neglected. The boundary condition at the
wellborne was given separately and independently for each layer, just as in
Boulton and Streltsova (1977a), thus avoiding the difficulties.
The mathematical model formulated in Liu and Chert (1987) not only describes
the fluid flow through a multiple-porosity, multiple-permeability medium, but also
the flow through a multiple-layered reservoir:

K 1 0 [raSP] 0so (23a)

rarrL -~-rJ+ A S p = B O--t-'

102 s3
f i r' ' T 9 i 9 r w ~ T i

101 u, : u~ :

10 0 ~. ,.I ~ _ _

, l . i

I ~ t/ /7/"

c a l c u lIa t e d by
r e s u l t s

10-3 ~ // the presented paper

/i // ~ r . = 0 . 0 2 , r. = I O O )
/; / ~ .... ~on,,.o? .or N e u , ~ o
10 .4 ~ , , . , .,' . . . . z ,. , . ~ . , . , , , . , tD
10 .4 10-3 10-2 10-1 10 0 101 10 2 10 3 10 4 10 5 10 6

Fig. 11. Comparison of the exact solution of Chen et al. with that of Neuman and Witherspoon
neglecting horizontal flow in the aquitard (after Chen et al., 1986).
176 Z.-X. CHEN

Os,~ = 0,
(23b)
Or F~

(23c)
k Or Jr~ r~

SD(rw, t)= SD~(t)I, SDw(t) is to be determined, (23d)

SD(r, 0) = 0, (23e)
where

0Lll a12 OQ N - "~IC l

0L21 --0/22 OL2N 4'2 c2 0

A= , B=/x

O/N1 OLN2 - - O~NN. ~bNcN

kl -sol 7 -1 0" I
k2 0 SD2 I 1 0
~

K= SD I= .

kN - SON -1. .0-

u 27rh
Olii = E OLij ~ O~ij ~--- Oiji ~ S D ~- "S~
j=l q/x

N is the total number of the porous systems (or layers).

A general, exact solution to this problem was obtained by Liu and Chert (1987)
through an orthogonal transformation. This transformation made the system of
equations of the resulting eigenvalue problem uncoupled, with each of them
being an independent Bessel equation of zero order9 A variety of problems of
practical importance can be obtained by specifying the elements of the matrix of
crossflow, A, and these problems, in principle, can be considered as solved.
Type curves which take the interlayer crossflow into account, were recom-
mended for analyzing buildup data from interference tests (Streltsova, 1984).
Prijambodo et al. (1985) emphasized that before the long-term performance of
the wen can be represented by the performance of an equivalent singteqayer
system, there is a very short early period where the reservoir behaves as a
commingled system9 There is also a transitional period when the pressure
response of the well depends on the contrast in permeability and on the degree of
communication between the layers9 Understanding the short-time behavior is
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 177

very important because the time span of virtually all pressure buildup tests
encompasses the duration during which a layered reservoir with crossflow may
not behave as if it were a single-layer system. The important influence of the
different skin factors of various layers on the well response was indicated in
Prijambodo et al. (1985).
Two papers (Gao, 1987; Ehlig-Economides, 1987) were devoted to the
determination of individual-layer parameters. All the methods to determine
properties of individual layers proposed in the two papers, as well as those
proposed in Kucuk et al. (1986) on the same topic, require the data from the
wellbore pressure and flow rate of each layer to be measured simultaneously. The
method of Gao (1987) was based on an approximate solution of the mathematical
model proposed by the author himself (Gao, 1984), while the method in Ehlig-
Economides (1987) relied on a Laplace space solution and its several limiting
forms.

15. Concluding Remarks

Great progress in the theory of fluid flow through naturally fractured reservoirs
has been made during the last three decades due to the continuum approach
proposed by Barenblatt and Zheltov (1960). The major contributions can be
summarized as in Table I.
The continuum approach is appropriate for cases where fractures are
developed somewhat uniformly and spread all over the reservoir, but this
approach ceases to be effective if fractures exist as enumerable, discrete net-

Table I. Historical review of solutions to problems of flow in naturally fractured reservoirs

Year Authors Contributions

1960 Barenblatt and Zheltov Mathematical model [Pseudosteady-state interporosity flow

(PIF)]
1960 Barenblatt et al. Line-source solution (D = o~ = 0, PIF)
1963 Warren and Root Asymptotic solution (D = 0, PIF); method to estimate k2, w
and A
1969 Kazemi et al. Line-source solution (D = 0, PIF)
1976 De Swaan Mathematical model with transient interporosity flow (TIF)
(D = 0); approximate analytical solution
1977 Boulton and Streltsova Line-source solution (D = 0, TIF)
1977 Jiang Exact solution (D = 0, PIF)
1979 Mavor and Cinco-Ley Introduction of Stehfest's algorithm; solution including CD
and s (D = 0, PIF)
1980 Chen and Jiang Exact solution (PIF)
1980 Bourdet and Gringarten Type curves using independent variables tD/CD, C D e 2s, r
and h e 2s (D = 0 ; PIF, TIF)
1983 Streltsova Clarification of existence of semilog straight-line segment
Serra et al. during transition period (D = 0, TIF)
Cinco-Ley and Samaniego
1987 Liu et al. Exact solution including CD and sl (i = 1, 2) (PIF)
178 z.-x. CHEN
works. In this case there is ample scope for ability of the discrete approach, and it
is certainly worth developing. Some interesting achievements obtained recently
at the Lawrence Berkeley Laboratory (Long et al., 1982; Long and Billaux,
1987) and progress in the stochastic modeling of dispersion in fractured media
made by Schwartz and Smith (1987) present promising means of solving this very
difficult problem.
To solve the problem of fluid flow through double-porosity, double-per-
meability media, the Laplace transformation is commonly used. The solution in
Laplace space is usually easy to obtain, but its inversion to real space is
frequently difficult. So far only a limited number of exact solutions have been
obtained. After the introduction of the Stehfest algorithm into the reservoir flow
investigations in 1978, the numerical inversion became commonplace. In con-
trast, the method of separation of variables developed and used by Chen, his
colleagues and students demonstrates its impressive strength in obtaining exact
solutions in such cases. It can be expected that new exact solutions for other
complex flow problems might be obtained by the use of this method. However,
the method of Laplace transformation has the advantage of obtaining simple
asymptotic solutions which are quite convenient and important for practical use.
With the increasing complexity of the solutions associated with the increasing
difficulty of the problem, the possibility of using conventional methods of analysis
to estimate reservoir parameters is reduced. Even the dimensionless type-curve
matching method loses its usefulness because there are too many parameters to
be treated. Automated matching with least-squares optimization may be accepted
as an appropriate alternative. But for a computerized search for the minimum of
a function, very often it may not be possible to locate the minimum because the
convergence of search is dependent on the nature of the function and an initial
guess. A modified Horner analysis method proposed recently by Chen (1988)
overcomes the nonconvergence problem associated with a computerized search,
and this method appears promising.
Even though progress has been made in the theory of fluid flow through
double-porosity, double-permeability systems, there still remains a number of
important unsolved problems. The mechanism of interporosity flow needs further
investigation. No exact solution has been obtained for double-porosity, double-
permeability media with a transient interporosity flow regime. No transient well
test analysis method to estimate individual properties of double-porosity, double
permeability systems is available. These are matters that should be considered.
Besides the many unsolved linear problems such as those listed above, there
are also many nonlinear problems. When the linear Darcy's law is no longer valid
due to high flux through the fractures with wide aperture, and/or when the fluid is
highly compressible, e.g. gas, the problem of single-phase flow through double-
porosity, double-permeability media must be a nonlinear one. Solving such
nonlinear problems will be an important task for future research.
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 179

Acknowledgments
This work was done during my stay as a visiting scientist at the Division of
Petroleum Engineering and Applied Geophysics, the Norwegian Institute of
Technology, Funds were provided by Statoil, Norsk Hydro, and the Royal
Norwegian Council for Scientific and Industrial Research (NTNF). I would like
to thank Curtis H. Whitson for his valuable suggestions and help.

References
Abdassah, D. and Ershaghi, I., 1986, Triple-porosity systems for representing naturally fractured
reservoirs, Soc. Pet. Eng. J. April, 113-127; Trans. A I M E 281.
Aguilera, R., 1987a, Well test analysis of naturally fractured reservoirs, SPEFE, Sept., 239-252.
Aguilera, R., 1987b, Multiple-rate analysis for pressure-buildup tests in reservoirs with tectonic,
regional, and contractional natural fractures, SPEFE, Sept., 253-260.
Avagan, E. A., 1967a, Some approximate solutions to problems of flow in fracture-pore medium, Izv.
Akad. Nauk, USSR, Mekh. Zidk. i Gasa, No. 4, 108-113 (in Russian).
Avagan, E. A., 1967b, Axisymmetrical problems of unsteady state flow of fluids in naturally fractured
reservoirs, Trudi VNII, No. 50, Nedra, Moscow, pp. 11-22 (in Russian).
Baker, W. J., 1955, Flow in fissured formations, Proc. Fourth World Petroleum Congress, Section II,
379-393.
Ban, A., 1961, Determination of time lag of build-up of pressure in fractured reservoirs, Izv. Akad.
Nauk, USSR, Otdel Tech. Nauk, Mekh. i Mash. No. 4, 38-42 (in Russian).
Barenblatt, G. E., Barisov, U. P., Kamenetsky, C. G., and Krilov, A. P., 1957. Determination of
parameters of oil reservoirs according to buildup data in shut-in wells, Izv. Akad. Nauk, USSR,
Otdel Tech. Nauk, No. 11 (in Russian).
Barenblatt, G. I. and Zheltov, Y. P., 1960, On fundamental equations of flow of homogeneous
liquids in naturally fractured rocks, Dokl. Akad. Nauk, USSR 132, No. 3, 545-548 (in Russian).
Barenblatt, G. I., Zheltov, Y. P., and Kochina, I. N., 1960, Basic concepts in the theory of seepage of
homogeneous liquids in fractured rocks, Prik. Matem. i Mekh. 24, No. 5,852-864 (in Russian).
Barenblatt, G. I., 1963, On certain boundary-value problems for the equations of seepage of a liquid
in fractured rocks, Prik. Matem. i Mekh. 27, No. 2, 348-350 (in Russian).
Bondarev, E. A. and Nicolaevsky, V. N., 1966, Formulation of problems in theory of seepage of
homogeneous liquid in naturally fractured porous media, N T S po Dabiche Neft, VNII, No. 30,
Gostoptechizdat, Moscow, pp. 29-33 (in Russian).
Boulton, N. S. and Streltsova, T. D., 1977a, Unsteady flow to a pumped well in a two-layered
water-bearing formation, J. Hydrology 35, 245-256.
Boulton, N. S. and Streltsova, T. D., 1977b, Unsteady flow to a pumped well in a fissured
water-bearing formation, J. Hydrology 35, 257-270.
Bourdet, D. and Gringarten, A. C., 1980, Determination of fissured volume and block size in
fractured reservoirs by type-curve analysis, Paper SPE 9293 presented at the 55th SPE Annual
Technical Conference and Exhibition, Dallas, TX, Sept. 21-24, 1980.
Bourdet, D., Ayoub, J. A., Whittle, T. M., Pirard, Y. M., and Kniazeff, V., 1983, Interpreting well
tests in fractured reservoirs, World Oil, Oct., 77-87.
Bourdet, D., 1985, Pressure behavior of layered reservoirs with crossflow, Paper SPE 13628
presented at the SPE 1985 California Regional Meeting, Bakersfield, California, March 27-29,
1985.
Braester, C., 1984, Influence of block size on the transition curve for a drawdown test in a naturally
fractured reservoir, Soc. Pet. Eng. J., Oct., 498-504; Trans. A I M E 277.
Charny, I. A., 1961, Flow of fluid in a formation with impermeable overburden and underlayer,
separated by a slightly permeable intercalation, Trudi M I N H i GP, No. 33, Gostoptechizdat,
Moscow, pp. 122-130 (in Russian).
Chen, C. C., Yeh, N., Raghavan, R., and Reynolds, A. C., 1984, Pressure response at observation
wells in fractured reservoirs, Soc. Pet. Eng. J., Dec., 628-638; Trans. A I M E 277.
180 Z.-X. CHEN

Chen, C. C., Serra, K., Reynolds, A. C., and Raghavan, R., 1985, Pressure transient analysis methods
for bounded naturally fractured reservoirs, Soc. Pet. Eng. J., June, 451-464; Trans. A I M E 279.
Chen, T., 1988, Pressure buildup analysis for well influenced by fractures and other heterogeneities,
J. Can. Pet. Tech., March-April, 61-71.
Chen, Z.-X. and Jiang, L.-S., 1980, Exact solution for the system of flow equations through a medium
with double-porosity, Scientia Sinica 28, No. 7,880-896.
Chen, Z.-X., 1982, Flow of homogeneous and heterogenous fluids through naturally fractured
reservoirs, Proc. International Meeting on Oilfield Development Techniques, Daqing, China, Sept.
1982, The Petroleum Industry Press, Beijing, 1985, pp. 332-354 (in Chinese).
Chen, Z.-X., 1983, Exact solution to the problem of flow of slightly compressible fluids in a bounded
confined 'fracture-pore' medium and its application to well testing, Applied Mathematics and
Mechanics (English Edition) 4, No. 3,435-477.
Chert, Z.-X., Pang, Z.-Y., Jiang, L.-S., and Liu, M.-X., 1986, Exact solution for the problem of
crossflow in a bounded two-aquifer system with an aquitard, Water Resour. Res. 22, No. 8,
1225-1236.
Chen, Z.-X. and You, J., 1987, The behavior of naturally fractured reservoirs including fluid flow in
matrix blocks, Transport in Porous Media 2, No. 2, 145-163.
Chen, Z.-X., 1988, Exact solutions for double-porosity, double-permeability systems including skin
effect for constant-pressure condition at the wellbore, Paper 6, Dr. Techn. thesis, Norwegian
Institute of Technology, Division of Petroleum Engineering and Applied Geophysics, 1988.
Cinco-Ley, H. and Samaniego, F. V., 1982, Pressure transient analysis for naturally fractured
reservoirs, Paper SPE 11026 presented at the 57th SPE Annual Fall Technical Conference and
Exhibition, New Orleans, LA, Sept. 26-29, 1982.
Closmann, P. J., 1975, An aquifer model for fissured reservoirs, Soc. Pet. Eng. J., Oct., 385-398;
Trans. A I M E 259.
Cobb, W. M., Ramey, H. J. Jr., and Miller, F. G., 1972, Well-test analysis for wells producing
commingled zones, J. Pet. Tech., Jan., 27-37; Trans. A I M E 257.
Crawford, G. E., Hagedorn, A. R., and Pierce, A. E., 1976, Analysis of pressure buildup tests in a
naturally fractured reservoir, J. Pet. Tech., Nov., 1295-1300.
Da Prat, G., Cinco-Ley, H., and Ramey, H. J. Jr., 1981, Decline curve analysis using type-curves for
two-porosity systems, Soc. Pet. Eng. J., June, 354-362; Trans. A I M E 271.
De Glee, G. J., 1930, Over Groundwaterstroomingen bij Wateronnrekking door Middel van Putten, T.
Waltman, Jr., Delft.
De Swaan, A., 1976, Analytic solutions for determining naturally fractured reservoir properties by
well testing, Soc. Pet. Eng. J., June, 117-122; Trans. A I M E 261.
Deruyck, B. G., Bourdet, D. P., Da Prat, G., and Ramey, H. J. Jr., 1982, Interpretation of
interference tests in reservoirs with double porosity behavior - theory and field examples, Paper
SPE 11025 presented at the 57th SPE Annual Technical Conference and Exhibition, New Orleans,
LA, Sept. 26-29, 1982.
Dontsov, K. M. and Boyrchuk, V. T., 1971, Effect of characteristics of fractured media on pressure
buildup behavior, Izv. Vyssh. Uchebn. Zaved., Neft i Gas, No. 1, 42-46 (in Russian).
Dougherty, D. E. and Babu, D. K., 1984, Flow to a partially penetrating well in a double-porosity
reservoir, Water Resour. Res. 20, No. 8, 1116-1122.
Duguit, J. O. and Lee, P. C. Y., 1977, Flow in fractured porous media, Water Resour. Res. 13, No. 3,
558-566.
Duvaut, G., 1961, Drainage des syst~mes htttrog~nes, Revue IFP 16, No. 10, 1164-1181.
Earlougher, R. C. Jr., Kersch, K. M., and Kunzman, W. J., 1974, Some characteristics of pressure
buildup behavior in bounded multiple-layered reservoirs without crossflow, J. Pet. Tech., Oct.,
1178-1186; Trans. A I M E 257.
Earlougher, R. C. Jr., 1977, Advances in Well Test Analysis, Monograph Series, Vol. 5, Society of
Petroleum Engineers of AIME, Dallas.
Ehlig-Economides, C. A., 1987, A new test for determination of individual layer properties in a
multilayered reservoir, SPEFE Sept., 261-283.
Ershaghi, I. and Aflaki, R., 1985, Problem in characterization of naturally fractured reservoirs from
well test data, Soc. Pet. Eng. J., June, 445-450; Trans. A I M E 279.
Gao, C.-T., 1984, Single-phase fluid flow in a stratified porous medium with crossflow, Soc. Pet. Eng.
J., Feb., 97-106; Trans. A I M E 277.
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 181

Gao, C.-T., 1987, Determination of parameters for individual layers in multilayer reservoirs by
transient well tests, SPEFE, March, 43-65.
Ge, J.-L. and Wu, Y.-S., 1982, The behavior of naturally fractured reservoirs and the technique for
well test analysis at constant pressure condition, Petroleum Exploration and Development 9, No. 3,
53-65 (in Chinese).
Gringarten, A. C., Burgess, T. M., Viturat, D., Pelissier, J., and Aubry, M., 1981, Evaluating fissured
formation geometry from well test data: a field example, Paper SPE 10182 presented at the 56th
SPE Annual Technical Conference and Exhibition, San Antonio, TX, Oct. 5-7, 1981.
Gringarten, A. C., 1984, Interpretation of tests in fissured and multilayered reservoir with double-
porosity behavior: theory and practice, J. Pet. Tech., April, 549-564; Trans. A I M E 277.
Hantush, M. S. and Jacob, C. E., 1955, Nonsteady radial flow in an infinite leaky aquifer, Eos Trans.
A G U 36, No. 1, 95-100.
Hantush, M. S., 1960, Modification of the theory of leaky aquifers, J. Geophys. Res. 65, No. 11,
3713-3725.
Horner, D. R., 1951, Pressure build-up in wells, Proc. Third World Pet. Cong., The Hague, Sec. II,
503-523.
Horner, D. R., Pressure behavior in a well producing from a number of different horizons, Shell Oil
Co. report.
Huitt, J. L., 1956, Fluid flow in simulated fractures, J. Am. Inst. Chem. Eng. 2, 259-264.
Hunt, B., 1985, Flow to a well in a multiaquifer system, Water Resour. Res. 21, No. 11, 1637-1642.
Hurst, W., 1934, Unsteady flow of fluids in oil reservoirs, Physics 5, January, 20-30.
Jacob, C. E., 1946, Radial flow in a leaky artesian aquifer, Eos. Trans. A G U 27, No. 2, 198-208.
Jacquard, P., 1960, Etude math6matique du drainage d'un r6servoir het6rog6ne, Rev. Inst. Franc.
Pdtrole 15, No. 10.
Jiang, J.-G., 1977, Exact solutions to problem of flow of fluid in naturally fractured reservoirs,
Mechanics, No. 4, 263-269 (in Chinese).
Katz, M. L. and Tek, M. R., 1962, A theoretical study of pressure distribution and fluid flux in
bounded stratified porous systems with crossflow, Soc. Pet. Eng. J., March, 68-82; Trans. A I M E
225.
Kazemi, H., 1969, Pressure transient analysis of naturally fractured reservoirs with uniform fracture
distribution, Soc. Pet. Eng. J., Dec., 451-462; Trans. A I M E 246.
Kazemi, H., Seth, M. S., and Thomas, G. W., 1969, The interpretation of interference tests in
naturally fractured reservoirs with uniform fracture distribution, Soc. Pet. Eng. J., Dec., 463-472;
Trans. A I M E 246.
Kazemi, H., 1970, Pressure buildup in limit testing of stratified systems, J. Pet. Tech., April,
503-511; Trans. A I M E 249.
Kucuk, F., Karakas, M., and Ayestaran, L., 1986, Well testing and techniques for layered reservoirs,
SPEFE, Aug., 342-354.
Kutliarov, V. S., 1967, Determination of parameters of naturally fractured reservoirs from data of
unsteady state flow of liquid towards wells, Trudi VNII, No. 50, Gostoptechizdat, Moscow, pp.
109-117 (in Russian).
Larsen, L., 1981, Wells producing commingled zones with unequal initial pressures and reservoir
properties, Paper SPE 10325 presented at the 56th SPE Annual Technical Conference and
Exhibition, San Antonio, TX, Oct. 5-7, 1981.
Larsen, L., 1982, Determination of skin factors and flow capacities of individual layers in two-layered
reservoirs, Paper SPE 11138 presented at the 57th SPE Annual Technical Conference and
Exhibition, New Orleans, TX, Sept. 26-29, 1982.
Lai, C. H., Bodvarsson, G. S., Tsang, C. F., and Witherspoon, P. A., 1983, A new model for well test
data analysis for naturally fractured reservoirs, Paper SPE 11688 presented at the 1983 SPE
California Regional Meeting, Ventura, CA, March 23-25, 1983.
Lang, Z.-X. and Li, Y., 1981, Identification of reservoir parameters in double-porosity systems, J.
East China Petroleum Inst., No. 3, 6-15 (in Chinese).
Lefkovits, H. C., Hazebroek, P., Allen, E. E., and Matthews, C. S., 1961, A study of the behavior of
bounded reservoirs composed of stratified layers, Soc. Pet. Eng. J., March, 43-58; Trans. A I M E
222.
Liu, C. Q., 1981, Exact solution for the compressible flow equation through a medium with triple
porosity, Applied Mathematics and Mechanics 2, No. 4, 457-461.
182 z . - x . CHEN

Liu, C.-Q. and An, W.-T., 1982, Numerical simulation of flow of a slightly compressible fluid through
a medium with multiple porosity, Acta Mechanica Sinica 14, No. 3, 236-243 (in Chinese).
Liu, M.-X. and Chen, Z.-X., 1987, A general solution for flow of fluid through multiple-porosity or
multiple-layered media, unpublished.
Liu, X.-N., Chen, Z.-X., and Jiang, L.-S., 1987, Exact solution of double-porosity, double-per-
meability systems including wellbore storage and skin effect, Paper SPE 16849 presented at the
62nd SPE Annual Conference and Exhibition, Dallas, TX, Sept. 27-30, 1987.
Lomize, G. M., 1951, Seepage Through Fractured Rocks, Gosenergoizdat, Moscow (in Russian).
Long, J. C. S., Remer, J. S., Wilson, C. R., and Witherspoon, P. A., 1982, Porous media equivalents
for networks of discontinuous fractures, Water Resour. Res. 18, No. 3, 645-658.
Long, J. C. S. and Billaux, D. M., 1987, From field data to fracture network modeling: an example
incorporating spatial structure, Water Resour. Res. 23, No. 7, 1201-1216.
Luan, Z.-A., 1981a, Problem of fluid flow through a double-porosity medium including the effects of
afterflow and skin, J. South-Western Petroleum College, No. 1, 50-58 (in Chinese).
Luan, Z.-A., 1981b, Analytical solution for transient flow of non-Newtonian fluids in naturally
fractured reservoirs, Acta Petrolei Sinica 2, No. 4, 75-79 (in Chinese).
Maksimov, V. A., 1960, The influence of nonhomogeneities on the determination of reservoir
parameters based on data of unsteady-state fluid influx into wells. A case of a two-layer formation,
Izv. Akad. Nauk, USSR, Otdel Tech. Nauk, Mekh. i. Mash., No. 3 (in Russian).
Matthews, C. S. and Russell, D. G., 1967, Pressure Buildup and Flow Tests in Wells, Monograph
Series, Vol. 1, Society of Petroleum Engineers of AIME, Dallas.
Mavor, M. J. and Cinco-Ley, H., 1979, Transient pressure behavior of naturally fractured reservoirs,
Paper SPE 7977 presented at the 54th SPE California Regional Meeting, Ventura, CA, April
18-20, 1979.
McGuinness, M. J., 1986, Pressure transmission in a bounded randomly fractured reservoir of
single-phase fluid, Transport in Porous Media 1, 371-397.
McNabb, A., 1978, Factorizable 'fruitcake' boundary value problem, Supplement to Newsletter 12,
The New Zealand Mathematical Society, 9-18.
Miller, C. C., Dyes, A. B., and Hutchinson, C. A. Jr., 1950, The estimation of permeability and
reservoir pressure from bottom hole pressure build-up characteristics, Trans. A I M E 189, 91-104.
Moench, A. F., 1984, Double-porosity models for a fissured groundwater reservoir with fracture skin,
Water Resour. Res. 20 No. 7, 831-846.
Muskat, M., 1934, The flow of compressible fluids through porous media and some problems in heat
conduction, Physics 5, March, 71-94.
Muskat, M., 1937a, The Flow of Homogeneous Fluids Through Porous Media, McGraw-Hill, New
York.
Muskat, M., 1937b, Use of data on build-up hole pressures, Trans. A I M E 123, 44-48.
Muskat, M., 1949, Physical Principles of Oil Production, McGraw-Hill, New York.
Najurieta, H. L., 1980, A theory for pressure transient analysis in naturally fractured reservoirs, 3.
Pet. Tech., July, 1241-1250; Trans. A I M E 269.
Neuman, S. P. and Witherspoon, P. A., 1969, Theory of flow in a confined two-aquifer system, Water
Resour. Res. 5, No. 4, 803-816.
Odeh, A. S., 1965, Unsteady state behavior of naturally fractured reservoirs, Soc. Pet. Eng. J.,
March, 60-64; Trans. A I M E 234.
Ozkan, E., Ohaeri, U., and Raghavan, R., 1987, Unsteady flow to a well produced at a constant
pressure in a fractured reservoir, SPEFE, June, 186-200.
Papadopulos, I. S., 1966, Nonsteady flow to multiaquifer wells, 3". Geophys. Res. 71, No. 20,
4791 4797.
Parsons, R. W., 1966, Permeability of idealized fractured rock, Soc. Pet. Eng. J, June, 126-136;
Trans. A I M E 246.
Pendergrass, J. D. and Berry, V. J. Jr., 1962, Pressure transient performance of a multilayered
reservoir with crossflow, Soc. Pet. Eng. J., Dec., 347-354; Trans. A I M E 225.
Pirson, R. S. and Pirson, S. J., 1961, An extension of lhe Pollard analysis method of well pressure
build-up and drawdown tests, Paper SPE 101 presented at the 36th SPE Annual Meeting, Dallas,
TX, Oct. 8-11, 1961.
Pollard, P., 1959, Evaluation of acid treatments from pressure build-up analysis, Trans. A I M E 216,
38-43.
TRANSIENT FLOW OF SLIGHTLY COMPRESSIBLE FLUIDS 183
Pottier, J., 1961, Mod61es h6t6rog6nes avec communications, Rev. Inst. Fran. P~trole 16, No. 10,
1127-1146.
Prijambodo, R., Raghavan, R., and Reynolds, A. C., 1985, Well test analysis for wells producing
layered reservoirs with crossflow, Soc. Pet. Eng. J., June, 380-396; Trans. A I M E 279.
Qi, Y.-F. and Zhang, X., 1986, Parameter estimation in a naturally fractured reservoir by the method
of cybernetics, Paper SPE 14864 presented at the SPE 1986 International Meeting on Petroleum
Engineering, Beijing, China, March 17-20, 1986.
Raghavan, R., Topaloglu, H. N., Cobb, W. M., and Ramey, H. J. Jr., 1974, Well-test analysis for
well producing from two commingled zones of unequal thickness, J. Pet. Tech., Sept., 1035-1041;
Trans. A I M E 257.
Reynolds, A. C., Chang, W. L., Yeh, N., and Raghavan, R., 1985, Wellbore pressure response in
naturally fractured reservoirs, J. Pet. Tech., May, 908-920; Trans. A I M E 279.
Romm, E. S., 1966, Flow Characteristics in Fractured Rocks, Nedra, Moscow (in Russian).
Rubinstein, L. I., 1948, Process of conduction of heat in heterogeneous media, Izv. Akad. Nauk,
USSR, Ser. Geogr. i Geofiz. 12, No. 1, 27-45 (in Russian).
Russell, D. G. and Prats, M., 1962a, Performance of layered reservoirs with crossflow - single-
compressible-fluid case, Soc. Pet. Eng. J., March, 53-67; Trans. A I M E 225.
Russell, D. G. and Prats, M., 1962b, The practical aspects of interlayer crossflow, J. Pet. Tech., June,
589-594.
Schwartz, F. W., Smith, L., and Crowe, A. S., 1983, A stochastic analysis of macroscopic dispersion
in fractured media, Water Resour. Res. 19, No. 5, 1253-1265.
Schwartz, F. W. and Smith, L., 1987, An overview of the stochastic modeling of dispersion in
fractured media, Advances in Transport Phenomena in Porous Media, NATO ASI Series, Martinus
Nijhoff, Dordrecht, Boston, Lancaster, pp. 727-750.
Serra, K., Reynolds, A. C., and Raghavan, R., 1983, New pressure transient analysis methods for
naturally fractured reservoirs, Soc. Pet. Eng. J., Dec., 2271-2283; Trans. A I M E 275.
Shalimov, B. V., 1966, On interference of wells in double-porosity reservoirs, Izv. Akad. Nauk,
USSR, Mekh. Zidk. i Gasa No. 4, 143-147 (in Russian).
Shapiro, A. M., 1987, Transport equations for fractured porous media, Advances in Transport
Phenomena in Porous Media, NATO ASI Series, Martinus Nijhoff, Dordrecht, Boston, Lancaster,
pp. 405-471.
Snow, D. T., 1969, Anisotropic permeability of fractured media, Water Resour. Res. 5, No. 6,
1273.
Stehfest, H., 1970, Algorithm 368, numerical inversion of Laplace transforms, Commun. AC1 13,
No. 1, 47-49.
Streltsova, T. D., 1983, Well pressure behavior of a naturally fractured reservoir, Soc. Pet. Eng. J.,
Oct., 769-780; Trans. A I M E 275.
Streltsova, T. D., 1984, Buildup analysis for interference tests in stratified formations, J. Pet. Tech.,
Feb., 301-310; Trans. A I M E 277.
Sun, S.-H. and Bai, D.-H., 1978, A problem of identification in the theory of fluid flow through
naturally fractured reservoirs, J. Sichuan Univ. (Natural Sci. Edn), No. 1 (in Chinese).
Tariq, S. M. and Ramey, H. J. Jr., 1978, Drawdown behavior of a well with storage and skin effect
communicating with layers of different radii and other characteristics, Paper SPE 7453 presented at
the 53rd SPE Annual Technical Conference and Exhibition, Houston, TX, Oct. 1-3, 1978.
Tempelaar-Lietz, W., 1961, Effect of oil production rate on performance of wells producing from
more than one horizon, Soc. Pet. Eng. Y., March, 26-31; Trans. A I M E 222.
Theis, C. V., 1935, The relation between the lowering of the piezometric surface and the rate and
duration of discharge of a well using ground-water storage, Trans. Amer. Geophys. Union 16,
519-524.
Uldrich, D. O. and Ershaghi, I., 1979, A method for estimating the interporosity flow parameter in
naturally fractured reservoirs, Soc. Pet. Eng. J., Oct., 324-332; Trans. A I M E 267.
Van Golf-Racht, T. D., 1982, Fundamentals of Fractured Reservoir Engineering, Developments in
Petroleum Science, 12, Elsevier Scientific, Amsterdam, Oxford, New York, p. 365.
Warren, J. R. and Root, P. J., 1963, The behavior of naturally fractured reservoirs, Soc. Pet. Eng. J.,
Sept., 245-255; Trans. A I M E 228.
Wilson, C. R. and Witherspoon, P. A., 1974, Steady state flow in rigid networks of fractures, Water
Resour. Res. 10, No. 2, 328-335.
184 z . - x . CHEN

Wu, Y.-S. and Ge, J.-L., 1981, Transient flow in a composite naturally fractured reservoir, Petroleum
Exploration and Development 9, No. 4, 58-67 (in Chinese).
Wu, Y.-S. and Ge, J.-L., 1983, The transient flow in naturally fractured reservoirs with three-porosity
systems, Acta Mechanica Sinica 26, No. 1, 81-85 (in Chinese).
You, J. and Chen, Z.-X., 1986, An improvement in the determination of interporosity flow parameter
in a double-porosity medium, Petroleum Exploration and Development 13, No. 3, 56-62 (in
Chinese).
Zhang, H.-Q., 1982, An analytical solution of gas flow towards a partially penetrated well in double
porous media, Acta Petrolei Sinica 3, No. 2, 51-62 (in Chinese).
Zhu, Y.-D., Zhang, J.-X., Le, G.-M. and Shi, J.-H., 1981, Study of the unsteady state flow of slightly
compressible fluid in naturally fractured reservoirs and its application, Aeta Petrolea Sinica 2, No.
3, 63-72 (in Chinese).