Flowing Bottomhole Pressure Calculation for a Pumped Well under Multiphase Flow.

Authors: Bikbulatov S., Khasanov M., Zagurenko A.

Summary
The ability to monitore bottomhole flowing pressure in pumping oil wells provides important information
regarding both reservoir and artificial lift performance. Converting surface pressure measurements to bottomhole is
currently accomplished by locating the fluid level in the annulus using a sonic device and then applying a correlation
to estimate the density of the gas-cut liquid column above the perforations.
This work proposes a flowing bottomhole pressure calculation procedure from fluid level measurements.
The model is developed from experimental work and from theoretical arguments.
The calculation procedure developed allows to calculate BHP without shutting the well, which is common
for fluid level casinghead pressure measurements. Also this method allows to take into account real geometry of the
well.
The comparison of the calculated and measured pump intake pressure shows good accuracy of a technique.
It allows to draw a conclusion for an opportunity to use this method in practice.

Significant progress has been made in the

Introduction acoustic device’s ability to monitor the movement of
the liquid column as function of time. Estimation of
change in liquid-column density, however, still
fraught with uncertainties.
"Is the well producing all the fluid that it is The fluid distribution in the annulus is a
capable of producing without problems?" function of the producing conditions of the particular
A. L. Podio well. The situation is generally found in the field: the
liquid level is above the formation and casinghead
The producing-rate efficiency of a well can gas is produced. This condition results in a gaseous
be determined with the curve of inflow performance annular liquid column. At stabilized producing
relationship1, which requires knowledge of the BHP. conditions, the oil in the casing annulus becomes
The ability to monitor bottomhole pressures saturated with the gas that is continuously flowing to
provides many advantages for reservoir management. the surface.
Pumping wells completed without a packer provide a The flow of gas through a static liquid
special opportunity for this low cost and reliable column creates a special type of multiphase flow,
bottomhole pressure surveillance. For these pumping termed Zero Net Liquid Flow (ZNLF) (fig.1). In this
wells, wellhead pressure data is converted to case, liquid is present in
bottomhole pressure by use of flow models and an the wellbore but does
acoustic (sonic) device to locate the gas-liquid not flow out the
interface. Fig. 1 shows the schematic of this type of tubing/casing annulus
pumping oil well. The well is completed in a with the gas phase. The
conventional fashion, without a packer. The pump gas phase simply passes
can be a sucker rod pump, PC pump or ESP. The pro- through the column of
duced fluids are pumped from the well through the liquid. A number of
tubing string, while produced and solution gas travels methods have been
up the tubing/casing annulus and is produced as developed to predict the
casinghead gas at the surface. Acoustic devises are liquid holdup of this
used to determine the depth to the gas-liquid inter- gas-cut liquid column.
face. Once it has been located, bottomhole pressure is The well-known
estimated through use of flow models to calculate the Gilbert2 chart correlates
pressure drop through the gas phase above the a Liquid Correction
interface and oil-water-gas mixture that exists below Factor (LCF). The LCF corrects the liquid density to
the interface. account for the effect of the gas phase. In modern
From the knowledge of the lengths of gas terms, the LCF is related to liquid holdup, HL, by:
and liquid columns, BHP can be estimated by adding (LCF )ρ L = ρ m = ρ L H L + ρ G (1 − H L ) (1)
the pressures exerted by these columns to the
In 1977, Godbey & Dimon3 presented a correlation of
casinghead pressure. Although simple in concept, this
vsg with the gas void faction, 1 – HL. The most widely
indirect calculation presents potential problems in
a use method was proposed by Podio et al.4 in 1980.
two areas: resolution of the acoustic device
Hasan et al.5 presented a model for prediction of the
measuring the gas/liquid interface and estimation of
LCF which allows a variety of fluid properties to be
the gas-entrained-liquid-column density.
considered. In 1994, Kabir & Hasan6 presented a
1
comprehensive review of the current methods and vSg (5)
discussed the accuracy of the various available HL = 1−
Co vm + vs
methods. Research in two areas - mass transfer in
where C0 and vs are given by Eq. 3 and 4,
gas/liquid systems and two-phase flow - has
respectively.
produced a wealth of information for predicting gas
Slug Flow. The drift-flux model of Eq. 5
void fractions in stagnant liquid columns. None of
also was applied in slug flow, but with different
these correlations, however, account for the effect of
values for C0 and vs given by C0 = 1.2 and
casing and tubing diameters.
The purpose of this project is to explore the vs = 0.35 gd L
(ρ − ρ g ) sinθ (1 + cosθ )1.2 (6)
relevant literature and develop a model for flowing ρL
bottomhole pressure calculation from fluid level
measurements using a practical range of flow 1.2 Flow model for tubing/casing annulus.
conditions encountered in a pumping-well annulus. Bubbly Flow. For bubbly flow in an oilwell
annulus, Hasan and Kabir propose that the void
1 Flow model fraction may be estimated by an expression of the
Various flow models are used for casing and form given by:
tubing/casing annulus description because of v gs
fg = (7)
different geometries. Works of Hasan and Kabir20,21,22 Av gs + B0 v∞
was taken as a basis for this flow models.
Hasan and Kabir also developed a For gas bubbling through an annulus, the
mechanistic model to predict pressure gradients in values of the parameters A and B0 are likely to
wellbores. To model flow-pattern transitions, Hasan depend on the ID and OD. A to be linearly dependent
and Kabir adapted an approach very similar to that of on the ID-to-OD ratio:
Taitel et al.23 Hasan and Kabir identified the same A = A0 + A1 (d t / d c ) (8)
four flow patterns: bubble flow, slug flow, chum Eq. 8, of course, needs to be verified by
flow, and annular flow. experimental work.
With the values of A0, A1 and В0
1.1 Flow model for casing. determined from the experimental data, the proposed
Bubble/Slug Transition. Transition from model for gas void fraction during bubbly flow
bubble flow (the condition of small bubbles dispersed becomes and Using the Harmathy correlation with a
throughout the flow cross section) to slug flow (when constant of 1.50 to represent the effect of liquid
the bubble becomes large enough to fill the entire properties on v∝ and allowing for the liquid
cross section) requires agglomeration or coalescence. superficial velocity, vls, Eq. 7 becomes
Bubbles, other than very small ones, generally follow v gs
fg = +
a zigzag path when rising through a liquid. This (1.97 + 0.371d t / d c )(vgs + vls )
results in collisions among bubbles, with the
consequent bubble agglomeration and formation of v gs
+
1.5[gσ (ρ l − ρ g ) / ρ l2 ]
0.25
larger bubbles, which
sin θ
vSg = (Co vSL + vs ) (2) Slug Flow. The model for gas void fraction
4 − Co in slug flow becomes
C0 is the flow coefficient given by Eq. 3 v gs
fg = +
Co = 
1.2 if d < 0.12 m or if vSL > 0.02 m / s (3) (1.82 + 0.9d t / d c )(v gs + vls )
2.0 if d > 0.12 m and if vSL < 0.02 m / s v gs
+
increases with an increase in the gas flow rate. Hasan [0.30 + 0.22d t / d c ] g (d t − d c )(ρ l − ρ g ) / ρ l
et al.24 reported that a transition to slug flow is
expected at a void fraction of 0.25. By use of a drift- Bubbly/Slug-Flow Transition. Transition
flux concept, the transition then can be expressed by from bubbly to slug flow apparently occurs at a void
Eq. 2. fraction of about 0.2533,34. The validity of this
The terminal rise velocity of small bubbles transition criterion, however, remains to be
given by established for an annular geometry. In their work,
Taitel et al assumed that at the point of transition, slip
 gσ (ρ − ρ g ) 4
1

vs = 1.53 L L2 (4) (vg-vl) between the phases is equal to the terminal rise
ρL 
  velocity, v∝. Because of the existence of a velocity
The Harmathy25 expression, Eq. 4, is used and concentration profile, this assumption may be
for the slip or bubble-rise velocity. Transition to slug inappropriate. Instead, use of Eq. 7 or its equivalent
flow takes place at superficial gas velocities greater can be made for a nonstagnant liquid column:
than that given by Eq. 2. v gs (9)
fg =
A(v gs + vls ) + B0 v∞
Bubble Flow. In bubble flows the Because Eq. 9 is expected to apply to the
expression for holdup, HL is entire bubbly flow regime, we can equate fg=0.25 at
the point of transition:
2
 v gs the total gradient). Acceleration, however, is small
0.25 =
A(v gs + vls ) + B0 v∞ and can be neglected.
In general, the acceleration component can
or
be neglected during all but the annular-flow pattern.
1
v gs = (0.25 Avls + 0.25B0 v∞ ) Suggests that an accurate estimation of the liquid
1 − 0.25 A holdup is essential when computing the elevation
With the appropriate expressions for the component. This component accounts for most of the
parameters A and B, transition equation may be pressure drop occurring in the bubble- and slug-flow
written as patterns. Because of the different hydrodynamics in
1 each flow pattern, estimations of holdup, HL, in-situ
v gs = ×
1 − 0.25(1.97 + 0.371d t / d c ) mixture density, ρm, and friction factor, f, are made
{ [
× 0.25(1.97 + 0.371d t / d c )vls + 0.375 gσ (ρ l − ρ g ) / ρ l2 ]
0.25
} separately.

2 BHP Calculation
1.3 Pressure drop calculation.
In bubble flow the total pressure gradient Recalculation of pressure was made on each
can be written as the sum of the gravitation or step on depth; the temperature undertook average for
hydrostatic head (dp/dD)el, friction (dp/dD)f and each of three intervals: casing, tubing/casing annulus
acceleration (dp/dD)acc components. Thus, below fluid level and tubing/casing annulus above
 dp   dp   dp   dp  fluid level.
  =  +  +  =
 dD  l  dD  el  dD  f  dD  acc
2.1 Design procedure of calculation
fv 2 ρ dv
= ρ s g∆D sin θ + m s + vm ρ s m (10)
2d dD Process of calculation BHP occurs in two
where stages:
−0.2
ρ v d I. The first (direct) step in the solution of this
f C = 0.046 g Sg  (1 + 75λLC ) problem is the calculation of the annulus oil level
 µ 
 g  (dynamic level) depth Dl as a function of the bottom
qL
λLC = hole pressure Pwf (see, for example, curve A in fig.2)
qL − q g
qL + q g II. At the second (reverse) step we use the
vm = = vSL + vSg calculated curve Dl = f (Pwf ) to estimate the
A
To estimate total pressure gradient, Eq. 10 bottomhole pressure Pwfl corresponding to the
can be used with the mixture density calculated from measured dynamic level Dlm .
the liquid holdup estimated with Eq. 5. The friction
component can be computed by treating the mul- Dl A
tiphase mixture as a homogeneous fluid. Friction
factor, f, can be determined from the Moody diagram Dlm i-th calculated
for a Reynolds number defined as point
Dli
ρ v
N Re m = L m Pwf
µL Pwfe Pwfi
This was recommended by Govier and Aziz
because would not be too different from pm/ftm, and Figure 2. – The procedure of bottomhole pressure
the contribution of the friction component to the total estimation.
pressure gradient is very small.
In slug flow, as in bubble flow, the total Calculation procedure (direct step):
pressure gradient can be obtained by Eq. 10 by use of 1. Calculate the average temperatures T 1 , T 2 , T 3 :
Eqs. 5,6. The estimation of the friction component
presents some difficulty because some of the liquid 2. Calculate the bubblepoint pressure average
flows downward in a film around the Taylor bubble values in the casing ( T = T 1 ) and in the annulus
while most of the liquid flows upward in the liquid ( T = T 2 ):
slug. Wallis35 suggested that the wall shear stress
Pb (T ) = PbR 10 0.00164 (T −TR )
around the vapor bubble be ignored. With this
assumption, the friction pressure gradient becomes
 dp  f C vm2 ρ L H L
  = 3. Calculate the annulus gas pressure at the
 dD  f 2d measured dynamic level Dlm :
The product pLHL is very nearly equal to pm
for low-pressure systems, indicating the similarity in  0.0342 Dlγ g 
Pl = Pa exp  
evaluating the friction terms in slug and bubble flow.  T 3 z3 
The contribution of the friction component is no
longer negligible but is still small (typically 10% of

3
 ∆D = 5 m or 10 m). The iterative procedure is to
The value of z 3 is calculated by the trial- be stopped when Pj +1 ≤ Pl .
and-error procedure:
e) Interpolate between the last two values of P
z3 = 1 to obtain the depth Dli corresponding to
 0.0342 Dlγ g  pressure Pl :
Pl = P1 exp   Pj − Pl
 T 3 z3  Dli = D j + ∆D
Pj − Pj +1

(P1 + Pl ) where D j the depth where the iteration was

P3 = stopped.
2 Addition calculations are presented into
z 3 = f z (P3 ,T3 ) Appendix 2.
where f z (P, T ) = gas compressibility as a function of
pressure P and temperature T. 2.2 Data collection
For the purpose of this study, data were
4. Determine the bottom hole pressure values to be collected from technology regimes on Priobskoye
set during calculations: field. These files are formed every month and contain
(Pwf max − Pwf min ) the information on an operating mode of all
Pwfi = Pwf min + (i − 1) producing wells.
N −1
As candidates for BHP calculation wells
which have been chosen are equipped with pressure
where i = 1,2,..., N , N – the number of calculations
gauges on the pump intake.
needed to construct the function Dl = f (Pwf ) (usually The following data were selected for
N≈20). calculation:
For each value of Pwfi (i = 1,2,..., N ) the - depth (well, setting pump and fluid dynamic
level);
following procedure is used to calculate the - volumetric flow rate (oil and water);
corresponding value of dynamic level depth Dli . - physical and chemical properties of fluids (oil
and gas specific gravity, oil viscocity);
a) Select a depth increment ∆D : - casing and tubing diameter.
(Dw − D p )
∆D =
Np Conclusions

The value of ∆D should not be greater than The procedure of BHP calculation from
10 m. dynamic fluid level for pumping wells completed
without packer was developed during performance of
b) Starting with known pressure value, Pwfi , at the given project.
The main advantage of this method that it
depth Dw , calculate a pressure traverse by allows to calculate BHP without measuring dp/dt at
the iterative procedure surface as it is done in standard methods of BHP
~  dp 
Pj +1 = Pj −  
j

P =
(P
j
~
+ Pj +1 ) calculation. Hence allows to simplify measurements
1 on well and to calculate BHP at any moment of
 dD  l j+
2
2
production. Also this method allows to take into
1
 dp 
j+
2 account real geometry of the well.
Pj +1 = Pj − D  D j = Dw − j∆D For a practical substantiation of a technique
 dD  l bottomhole pressure and corresponding pressure on
where j = 0,1,2,..., N p , P0 = Pwfi , T = T1 . the pump intake has been calculated on 85 wells of
Priobskoye field. Results of comparison of the
c) Calculate at the pump intake ( j = N p ): Pup, calculated and measured pressures are shown in
figure 3. Comparison shows good accuracy of this
wgup , f gup , qoup , υ sup . technique (R2=0.76).
Then estimate:
a. separation coefficient E s ,
b. annulus gas mass rate wga .

d) Calculate a pressure traverse in the annulus

by the iterative procedure (b), where now
D j = D p − j∆D , j = 0,1,2,.... , Po = Pup ,
T = T2 , ∆D is to be selected properly (usually
4
 µor oil viscosity at reservoir conditions, mPa⋅s
15,0 Np the number of depth increments,

Measured pressure at pump intake, MPa

2
R = 0,7636 dimensionless
12,5 P pressure, MPa
Pa annulus gas pressure at the wellhead, MPa
10,0 Pb bubblepoint pressure at temperature T, MPa
Pb1 bubblepoint pressure at temperature T 1 ,
7,5 MPa
PbR bubblepoint pressure at reservoir
5,0 temperature, MPa
Pl gas pressure at the dynamic level, MPa
2,5 Ppc pseudocritical pressure, MPa
Psc pressure at standard conditions, MPa
0,0 Pup pump intake pressure, MPa,
0,0 2,5 5,0 7,5 10,0 12,5 15,0 Pwfe the estimated bottom hole pressure, MPa
Pwfi i-th value of the bottom hole pressure (set),
Calculated pressure at pump intake, MPa MPa
Figure 3. – Calculated and Measured pressure at ∆p pressure increment, MPa
intake comparison. P3 average pressure in the annulus gas cap,
MPa
Purpose of further work qt total flow rate, m3/day,
qga gas volumetric flow rate in the annulus,
It is necessary to specify possible directions m3/day
of further work for project improvement: qo, qw, qg flow rates (volumetric) of oil, water and gas
- use in calculations of real distribution of at given P and T, m3/day
temperature on a well; qosc oil (volumetric) flow rate at standard
- creation of more exact calculation model for conditions, m3/day
quantity of separated gas on the pump intake; qup oil flow rate at pump intake, m3/day
- specification of used correlations for PVT ρg gas density at P and T, kg/m3
properties ρL liquid density at given pressure P and
temperature T
Nomenclature ρo oil density at P and T, kg/m3
A flow area of conduit, m2, ρw water density at P and T, kg/m3
Aa annular cross section area, m2 Rsi initial solution GOR, m3/ m3, dimensionless
Ac casing inside cross section area, m2 Rs(P,T) solution gas/oil ratio (GOR) at given P and T,
A, As parameters, dimensionless m3/m3
As0, As1 parameters, dimensionless NRe m Reynolds number of mixture, dimensionless
A0, A1 parameters, dimensionless T temperature, °K
Bo oil FVF at pressure P and temperature T, T1 average temperature in casing below the
m3/m3 pump intake, °K
B0 parameter, dimensionless T2 average temperature of the oil-gas mixture in
B, Bs, parameters, m/s the annulus, °K
Co oil isothermal compressibility, MPa-1 T3 average temperature of the annulus gas cap,
C0 parameter, dimensionless
°K
C1 parameter, m/s
Tpc pseudocritical temperature, °K
d c, d t casing ID and tubing OD, m
deq equivalent diameter, m TR reservoir temperature, °K
Dli dynamic level corresponding to Pwi Tsc temperature at standard conditions, °K
(calculated), m Twh wellhead temperature, °K
Dlm the measured dynamic level, m vs gas slip velocity, the difference between the
Dp pump setting depth, m average gas and liquid velocities, m/s
Dw bottomhole depth, m vs up gas slip velocity at pump intake, m/s
∆D depth increment, m vg actual gas velocity, m/s
Es separation coefficient, dimensionless vsg superficial gas velocity, m/s
Fgc gradient correction factor, dimensionless vl actual liquid velocity, m/s
Fwo water-oil ratio, dimensionless vm mixture velocity, m/s
f friction factor, dimensionless vs1 superficial liquid velocity, m/s
fg gas fraction, dimensionless vsw gas slip velocity in the water-oil mixture,
fg up gas fraction at pump intake, dimensionless m/s
g gravitational acceleration constant, m/s2 v∞ terminal raise velocity for a single bubble in
gT temperature gradient, °K/m infinite medium, m/s
Hl liquid holdup, dimensionless
5
wg free gas mass flow rate in casing at pressure Upward Gas-Liquid Flow in Vertical Tubes,"
P of interest, kg/day, AlChEJ. (1980)26,345.
wga gas mass flow rate in the annulus, kg/day, 24. Hasan, A.R., Kabir, C.S., and Rahman, R.:
wup free gas mass flow rate in the casing at "Predicting Liquid Gradient in a Pumping-Well
suction conditions, kg/day Annulus," SPEPE February 1988) 113; Trans.,
γo oil specific gravity, dimensionless AIME, 285.
γg gas specific gravity, dimensionless 25. Harmathy.T.Z.: "Velocity of Large Drops and
z gas compressibility factor, dimensionless Bubbles in Media of Infinite or Restricted Extent,"
z3 average gas cap compressibility, AIChE J. (1960) 6,281.
dimensionless 26. Mashetfcar R.A.: "Bubble Columns" Brit.
σ surface tension, dynes/cm Chem. Eng. (Oct. 1970) 1297-1304.
θ inclination angle of well 27. Haug H.F.: "Stability of Sieve Trays with High
Overflow Weirs." Chem. Eng. Sci. (1976) 31, 295.
SI Metric Conversion Factors 28. J. and Kastanek F.: "Gas Holdup in Uniformly
cP x 1.0*E-03 = Pa·s Aerated Bubble Column Reactors." Chem, Eng.
dyne x 1.0*E-02 = mN Commun. (1979) 3, Nos. 4 and 5, 413-29.
ft x 3.048* E-01 = m 29. Zuber H. and Hench J.: "Steady State and
o
F(°F-32)/1.8=°C Transient Void-Fraction of Bubbling Systems and
gal x 3.785 412 E-03 = m3 Their Operating Limits," Report No.62 GL 100,
in. x 2.54 = cm General Electric Co., Schenectady, NY (1962).
lbm x 4.535924E-01=kg 30. Wallis, G.B.: One Dimensional Two-Phase
psi x 6.895 = kPa Flow, McGraw Hill Book Co., New York City
o
API 141.5/(131.5 + °API) = g/cm3 (1969) 93, 248, 255.
ft x 3.048*E-01 = m 31. Nicklin D.J., Wilkes, J.O., and Davidson, J.F.:
psi x 6.895 = kPa "Two Phase Flow in Vertical Tubes," Trans., Inst.
Chem. Engrs. (Feb. 1962) 61-68.
References 32. Griffith, P.: "The Prediction of Low-Quality
. Vogel J.V.: "Inflow Performance Relationships Boiling Voids," Trans., ASME J. Heal Transfer
for Solution-Gas Drive Wells," JPTQw. 1968) 83- (Aug. 1964) 86, 327-33.
87; Trans., AIME, 243. 33. Taitel Y. and Dukler A.E.: "Flow Regime
2. Gipson, S.W. and Swaim, H.W.: "Designed Transition for Vertical Upward Gas-Liquid Flow:
Beam Pumping,' ' Proc. , 19th Annual A Preliminary Approach Through Physical
Southwestern Petroleum Short Course, Lubbock; Modeling," paper presented at the 1977 AIChE
TX (April 1972)95. Annual Meeting, New York City.
3. Godbey, J.K. and Dimon, C.A.: "The Automatic 34. Taitel Y., Bamea D., and Dukler A.E.:
Liquid Level Monitor for Pumping Wells," JPT "Modelling Flow Pattern Transition for Steady
(Aug. 1977) 1019-24. Upward Gas-Liquid Flow in Vertical Tubes."
4. Podio, A.L.: Tarrillion, ML, and Roberts, E.T.: AIChE J. (May 1980) 345-54, York City.
"Laboratory Work Improves Calculations," pit & 35. Steen, D.A and Wallis, G.B.: AEC report No.
Gas J. (Aug. 25, 1980) 137-46. NYO-31142-2 (1964).
5. Hasan, A.R. and Kabir, C.S.: "Determining
Bottomhole Pressures in Pumping Wells," SPEI
(Dec. 1985) 823-38.
6. Kabir, C.S. and A.R. Hasan: "Two-Phase Flow
Correlations as Applied to Pumping Well
Testing," ASME J. of Energy Resources Tech.,
Vol. 116, 121-127 (June 1994). October.
19. Ansari, A.M. et al.: "A Comprehensive
Mechanistic Model for Two-Phase Flow in
Wellbores." SPEPF (May 1994) 143; Trans.,
AIME, 297.
20. Hasan, A.R. and Kabir, C.S.: "A Study of
Multiphase Flow Behavior in Vertical Wells,"
SPEPE (May 1988) 263; Trans., AIME, 285.
21. Hasan, A.R. and Kabir, C.S.: "Predicting
Multiphase Flow Behavior in a Deviated Well,"
SPEPE (November 1988) 474.
22. Kabir, C.S. and Hasan, A.R.: "Performance of
a Two-Phase Gas/Liquid Model in Vertical
Wells," J. Pet. Sci. & Eng. (1990) 4,273.
23. Taitel. Y.M., Bamea, D.. and Dukler. A.E.:
"Modeling Flow Pattern Transitions for Steady
6