Site categories

Performance of the Network of Pipes in Gas Industry

Previous article Manual for Engineers to Calculate Amount of Gas in Reservoirs“Essential Calculations of the Gas Reservoir Performance” demonstrated that the ability of a gas reservoir to produce for a given set of reservoir conditions depends directly on the flowing bottom-hole pressure, pwf. This article will help engineers calculate gas pressure in piping system, flows in wells, pipelines and flowing temperature.

Picture Everything about Natural Gas in Modern Industry“Possible pressure losses in complete system” illustrates that pwf depends on the separator pressure and the configuration of the piping system; that is:

pwf = psep+Δpfl+Δpch+Δptub+Δpris       Eq. 1


  • psep = separator pressure;
  • Δpfl = pressure drop in the flowline;
  • Δpch = pressure drop in the surface choke;
  • Δptub = pressure drop in the well tubing, and;
  • Δpris = pressure drop in other restrictions, such as subsurface safety valves (SSSVS), valves and fittings, etc.

In order to determine the deliverability of the total well system, it is necessary to be able to calculate all of the pressure drops listed in Equation 1. This article will develop equations to make these calculations and demonstrate their application. The effects of liquids in the flow stream will be discussed, and the use of previously prepared pressure traverse curves for quick estimates will be demonstrated. Only steady-state, one-dimensional flow will be considered in this article.

Basic Flow Equation

The theoretical basis for most fluid flow equations is the general energy equation, an expression for the balance or conservation of energy between two points in a system. The energy equation is developed first, and using thermodynamic principles, is modified to a pressure gradient equation form.

The steady-state energy balance simply states that the energy of a fluid entering a control volume, plus any shaft work done on or by the fluid, plus any heat energy added to or taken from the fluid must equal the energy leaving the control volume. Picture 1 may be used to illustrate the control volume principle.

Control volume principle
Pic. 1 Flow system control volume

Considering a steady-state system, the energy balance may be written as:

U1+p1V1+mv122gc+m g h1gc+q+Ws=U2+p2V2+mv222gc+m g h2gc          Eq. 2


  • U

    = internal energy;

  • pV = energy of expansion or compression;
  • mv22gc

    = kinetic energy;

  • m ghgc

    = potential energy;

  • q

    = heat energy added to fluid, and;

  • Ws

    = work done on the fluid by the surroundings.

Dividing Equation 2 by m to obtain an energy per unit mass balance and writing in differential form gives:

dU=dpρ+vdvgc+ggcdh+dq+dWs=0.        Eq. 3

This form of the energy balance equation is difficult to pply because of the internal energy term, so it is usually converted to a mechanical energy balance using well-known thermodynamic relations. From thermodynamics:

dU=dhdpρ,          Eq. 4




dU=TdS+d pρdpρ          Eq. 5


  • h = enthalpy;
  • S = entropy, and;
  • T = temperature.

Substituting Equation 5 into Equation 3 and simplifying results in:

TdS+dpρ+vdvgc+ggcdh+dq+dWs=0            Eq. 6

For an irreversible process, the Clausius inequality states that:



TdS = -dq+dLw,


  • dLw = losses due to irreversibilities, such as friction.

Using this relationship and assuming no work is done on or by the fluid, Equation 6 becomes:

dpρ+vdvgc+ggcdh+dLw=0             Eq. 7

If we consider a pipe inclined at some angle θ to the horizontal, as in Picture 2, since dh = dL sinθ:

Scheme: flow geometry
Pic. 2 Flow geometry
dpρ+vdvgc+ggcdL sinθ+dLw=0

Multiplying the equation by



dpdL+ρvdvgcdL+ggcρ sinθ+ρdLwdL=0          Eq. 8

Equation 8 can be solved for pressure gradient, and if we consider a pressure drop as being positive in the direction of flow:

dpdL=ggcρ sinθ+ρvdvgcdL+dpdLf,             Eq. 9



is the pressure gradient due to viscous shear or friction losses.

In horizontal pipe flow the energy losses or pressure drops are caused by change in kinetic energy and friction losses only. Since most of the viscous shear occurs at the pipe wall, the ratio of wall shear stress (τw) to kinetic energy per unit volume (ρ v2/2 gc) reflects the relative importance of wall shear stress to the total losses. This ratio forms a dimensionless group and defines a friction factor.

fd=τwρv2/2gc=2τwgcρv2             Eq. 10

To evaluate the wall shear stress, a force balance between pressure forces and wall shear stress can be formed. Referring to Picture 3:

Scheme: force balance
Pic. 3 Force balance
p1p1dpdLdLπ d24=τw(πd)dL
τw=d4dpdLf                 Eq. 11

Substituting Equation 11 into Equation 10 and solving for the pressure gradient due to friction gives:


which is the well-known Fanning equation. In terms of a Darcy-Weisbach or Moody friction factor,


, and:

dpdLf=fρv22gcd,           Eq. 12

Laminar Single-Phase Flow

The friction factor for laminar flow can be determined analytically by combining Equation 12 with the Hagen-Poiseuille equation for laminar flow:




Equating the expressions for frictional pressure gradient gives:




The dimensionless group, NRe = ρvd/μ is the ratio of fluid momentum forces to viscous shear forces and is known as the Reynolds number. It is used as a parameter to distinguish between laminar and turbulent fluid flow. For engineering calculations, the dividing point between laminar and turbulent flow can be assumed to occur at a Reynolds number of 2 100 for flow in a circular pipe. Using units of lbm/ft3, ft/sec, ft and centipoise, the Reynolds number equation is:

NRe=1 488 ρvdμ

Turbulent Single-Phase Flow

The ability to predict flow behavior under turbulent flow conditions is a direct result of extensive experimental studies of velocity profiles and pressure gradients. These studies have shown that both velocity profile and pressure gradient are very sensitive to characteristics of the pipe wall. A logical approach to defining friction factors is to begin with the simplest case, i. e., the smooth-wall pipe, proceed to the partially rough wall, and finally to the fully rough wall. Only the most accurate empirical equations available for friction factors are presented here.

Smooth-Wall Pipe. For smooth-wall pipes, several equations have been developed, each valid over different ranges of Reynolds numbers. The most commonly used equation-since it is explicit in f and also covers a wide range of Reynolds numbers (3 000<NRe<3·106) – was presented by Drew, Koo, and McAdams in 1932.

f=0,0056+0,5NRe0,32           Eq. 13

An equation proposed by Blasius may be used for Reynolds numbers up to 100 000 for smooth pipes.

f=0,316NRe0,25           Eq. 14

Rough-Wall Pipe. The inside wall of a pipe is not normally smooth, and in turbulent flow, the roughness can have a definite effect on the friction factor, and thus the pressure gradient. Wall roughness is a function of the pipe material, the method of manufacture, and the environment to which it has been exposed.

From a microscopic sense, wall roughness is not uniform. Individual protrusions, indentations, etc. vary in height, width, length, shape, and distribution. The absolute roughness of a pipe, , is the mean protruding
height of relatively uniformly distributed and sized, tightly packed sand grains that would give the same pressure gradient behavior as the actual pipe.

Dimensional analysis suggests that the effect of roughness is not due to its absolute dimensions, but rather to its dimensions relative to the inside diameter of the pipe, ∈/d. In turbulent flow, the effect of wall roughness has been found to be dependent on both the relative roughness and on the Reynolds number. If the laminar sub-layer that exists within the boundary layer is thick enough, the behavior is similar to a smooth pipe. The sublayer thickness is directly related to the Reynolds number.

Nikuradse’s famous sand grain experiments formed the basis for friction factor data from rough pipes. His correlation for fully rough-wall pipe is still the best one available. The friction factor may be calculated explicitly from:

1f=1,742Log2d.             Eq. 15

The equation that is used as the basis for modern friction factor charts was proposed by Colebrook and White in 1939.

1f=1,742Log2d+18,7NRef             Eq. 16

The friction factor cannot be extracted readily from the Colebrook equation. By rearranging the equation as follows, a trial-and-error procedure may be used to solve the equation for friction factor.


Values of fg are estimated and then fc is calculated until fg and fc agree to an acceptable tolerance. Using the Drew, Koo, and McAdams equation as an initial guess is recommended. After each unsuccessful iteration, the calculated value becomes the assumed value for the next iteration. Also, if more than one pressure loss calculation is to be made as in the case of the iterative procedures discussed in later sections, then the “converged” value of the previous calculation should be used for the initial guess in the next calculation. Convergence using this method is rapid, normally taking only two or three iterations. The variation of single-phase friction factor with Reynolds number and relative roughness is shown graphically in Picture 4.

Variation of single-phase friction factor
Pic. 4 Friction factor for fully-developed flow in circular pipes

The Colebrook equation may be applied to flow problems in the smooth, transition, and fully rough zones of turbulent flow. For large values of Reynolds number, it degenerates to the Nikuradse equation.

An explicit friction factor equation was proposed by Jain and compared in accuracy to the Colebrook equation. Jain found that for a range of relative roughness between 10-6 and 10-2 and range of Reynolds number between 5·103 and 108 the errors were within ±1,0 % when compared with the Colebrook equation. The equation gives a maximum error of 3 % for Reynolds numbers as low as 2 000. The equation is:

1f=1,142Logd+21,25NRe0,9.            Eq. 17

Equation 17 is recommended for all calculations requiring a friction factor determination for turbulent flow. It is much easier to use than Equation 16 and, since the value of will usually not be known to any high degree of accuracy, will give satisfactory results.

The determination of the value to use for pipe wall roughness in the friction factor equations is sometimes difficult. It is important to emphasize that is not a property that is physically measured. Rather, it is the sand grain roughness that would result in the same friction factor. The only way this can be evaluated is by comparison of the behavior of a normal pipe with one that is sand-roughened. Moody has done this, and his results, given in Picture 5, are still the accepted values. These values should not be considered inviolate and could change significantly by such things as paraffin deposition, erosion, or corrosion. Thus, if measured pressure gradients are available, a friction factor and Reynolds number can be calculated, and an effective ∈/d obtained from the Moody diagram. This value of ∈/d obtained from the Moody diagram. This value of ∈/d should then be used for future predictions until updated again. If no information on roughness is available, a value of = 0,0006 ft is recommended for tubing and line pipe that has been in service for some time.

Example 1:

A liquid of specific gravity 0,82 and viscosity of 3 cp (.003 kg/m-sec) flows in a 4 in. (101,6 mm) diameter pipe at a velocity of 30 ft/sec (9,14 m/sec). The pipe material is new commercial steel. Caiculate the friction factor using both the Colebrook equation and the Jain equation.


From Picture 5, for commercial steel, ∈/d = 0,00045 Colebrook Solution: Use the Drew, Koo and McAdams equation for a first guess.

Roughness values of the commercial steel
Pic. 5 Relalive roughness values for pipes of common engineering materials
N = ρvd/μ = (820)(9,14)(.1016)/.003 = 253,824
fg = 0,015
fc = 0,0183

This value is not close enough to fg; therefore another trial is required using fg = 0,0183.

fc = 0,0182

A third trial using fg = 0,0182 gives fc = 0,0182. Jain Solution:

fc = 0,0183

Combining Equations 9 and 12, the pressure gradient equation, which is applicable to any fluid at any pipe inclination angle, becomes:

dpdL=ggcρ sin θ +fρv22gcd+ρvdvgcdL,             Eq. 18

where the friction factor, f, is a function of Reynolds number and pipe roughness. This relationship is shown in the Moody diagram (Pic. 4). The total pressure gradient can be considered to be composed of three distinct components; that is:

dpdL=dpdLel+dpdLf+dpdLacc,            Eq. 19


dpdLel=ggcρ sin θ

is the component due to potential energy or elevation change. It is also referred to as the hydrostatic component since it is the only component that would apply at conditions of no flow.


is the component due to friction losses.


is the component due to kinetic energy change or convective acceleration. Equation 18 applies for any fluid in steady-state, one-dimensional flow for which f, ρ, and v can be defined.

The elevation change or hydrostatic component is zero for horizontal flow only. It applies for compressible or incompressible, steady-state or transient flow in both vertical and inclined pipes. For downward flow the sine of the angle is negative, and the hydrostatic pressure increases in the direction of flow.

The friction loss component applies for any type of flow at any pipe angle. It always causes a drop of pressure in the direction of flow. In laminar flow the friction losses are linearly proportional to the fluid velocity. In turbulent flow the friction losses are proportional to vn, where 1,7≤n≤2.

The kinetic energy change or acceleration component is zero for constant area, incompressible flow. For any flow condition in which a velocity change occurs, such as compressible flow, a pressure drop will occur in the direction of the velocity increase.

Although single-phase flow has been studied extensively, it still requires an empirically determined friction factor for turbulent flow calculations. The dependence of this friction factor on pipe roughness, which must usually be estimated, makes the calculated pressure gradients subject to considerable error.

Equation 18 is a differential equation and must be integrated in order to apply it to calculate pressure drop as a function of flow rate or velocity and pipe diameter. It must be combined with a continuity equation and an equation of state to express velocity and density in terms of pressure. The following sections describe various assumptions made in integrating the equation for application.

If a computer is available, the equation can be integrated numerically by dividing the pipe into small increments and evaluating the gas or fluid properties at average pressure and temperature in the increment. If small enough increments are taken, the accuracy will be very good. A procedure for calculating a pressure traverse in a pipe using this method is outlined below, and a flow chart is presented in Picture 6.

Example of flow chart
Pic. 6 Flow chart for calculating a pressure traverse (incrementing on pressure)
  1. Starting with the known pressure, p1, at location L1 select a length increment, ΔL.
  2. Estimate a pressure increment, Δp, corresponding to length increment, ΔL.
  3. Calculate the average pressure and, for nonisothermal cases, the average temperature in the increment.
  4. From laboratory data or empirical correlations, determine the necessary fluid and PVT properties at conditions of average pressure and temperature (ρg, vg, μg).
  5. Calculate the pressure gradient, dp/dL, in the increment at average conditions of pressure, temperature, and pipe inclination, using Equation 18.
  6. Calculate the pressure increment corresponding to the selected length increment, Δp = ΔL (dp/dL).
  7. Compare the estimated and calculated values of Δp obtained in steps 2 and 6. If they are not sufficiently close, estimate a new pressure increment and return to step 3. Repeat steps 3 through 7 until the estimated and calculated values are sufficiently close.
  8. Set L = L1+ΣΔL and p = p1+ΣΔp.
  9. If ΣΔL is less than the total conduit length, return to step 2.

Using this procedure the length increments can be selected so that their sum is exactly equal to the total conduit length and interpolation is not required in the last step.

Flow in Wells

Several methods are available for calculating static and flowing pressure drop in gas wells. The most widely used method is that of Cullender and Smith. All of the methods begin with Equation 18, with modifications for flow geometry. In most cases the acceleration gradient is ignored. Since it is frequently necessary to calculate the static bottom-hole pressure in a gas well, this procedure will be presented first.

Static Bottom-Hole Pressure

Vor a vertical (θ = 90°, sin θ = 1), shut-in (v = 0) gas well, Equation 18 becomes:

dpdh=gρggc,              Eq. 20



Combining this with Equation (20):

dpp=gMdhgcZRT.           Eq. 21

Average Pressure and Temperature Method. If Z is evaluated at average pressure and temperature in the increment:

ptspwsdpp=gMgcR Z¯ T¯οHdh,

from which:

pws=ptsEXPgMHgcR Z¯ T¯,         Eq. 22

This equation holds for any consistent set of units. For conventional field units:

pws=ptsEXP0,01875γgH/T¯ Z¯           Eq. 23


  • pws = static or shut-in BHP, psia;
  • pts = static tubing pressure, psia;
  • γg = gas gravity (air = 1);
  • H = well depth, ft;
  • T¯

    = average temperature in the tubing, °R, and;

  • Z¯

    = gas compressibility factor evaluated at




Evaluation of


makes the calculation iterative, and the procedure outlined previously can be used.

Example 2:

Using the following data, calculate pws with Equation 23.

  • H = 10 000 ft;
  • Ts = 70 °F = 530 °R;
  • γg = 0,6;
  • pts = 4 000 psia;
  • Tf = 220 °F = 680 °R.


A good first guess for pws can be obtained from:

pws*=4 000(1+2,5·105(10 000))= 5 000 psia
T¯=Ts+Tf2=530+6802=605 °R
p¯=pts+pws*2=4 000+5 0002=4 500 psia
ppc = 709,6-58,7(.6) = 674,3 psia
Tpc = 170,5+307,3(.6) = 354,9 °R
ppr=p¯ppc=4 500/674,3=6,67

From “Best and Most Widely Methods to Perform Necessary Calculations with a Gas”“Compressibility factors for natural gases”, Z = 0,950.

pws = 4 000 EXP [(0,01875)(.6)(10 000)/(605)(Z)]
pws = 4 000 EXP [0,18595/Z] = 4 000 EXP [0,18595/0,950]
pws = 4 865 psia

This is not close enough to the estimated value of 5 000 psia. Set the calculated value of pws as the next estimated value and continue until convergence is reached.

pws*=4 865 psia   p¯=4 000+4 8652=4 433 psia
ppr = 4 433/674,3 = 6,57      Z = 0,947
pws = 4 000 EXP [0,18595/0,947] = 4 868 psia

This value is probably close enough to the previously estimated value of 4 865. However, one more iteration using 4 868 as the estimated value yields a calculated value of 4 868 also, which is therefore the correct bottomhole pressure.

The calculation could also be made by estimating an initial value for Z and comparing calculated and estimated values until convergence on Z is obtained.

Cullender and Smith Method. The method presented by Cullender and Smith takes into account the variation of temperature with depth and the variation of Z with pressure and temperature. From Equation 21:


The integral is written in short notation as:

ptspwsTZpdp=ptspwsI dp=0,01875γgH

Using a series expansion, the value of the integral is approximated by:

2  I dp=(pmspts)(Ims+Its)+(pwspms)(Iws+Ims),          Eq. 24


  • pms = pressure at mid-point of well, H/2;
  • Ims = I

    evaluated at pms,


  • Its = I evaluated at pts, Ts;
  • Iws = I evaluated at pws, Tf.

The calculation procedure consists of dividing the well into two equal segments of length, H/2, finding the pressure pms at H/2 and using this value to calculate pws. Its can be evaluated from known surface conditions; that is:


Example 3:

Work Example 2 using the Cullender and Smith method.


Calculate Its:

At T = 70, p = 4 000, Z = 0,86:

Its=TZp=530(.86)4 000=0,1140

Estimate pms:

pms*=4 000(1+2,5·105(5 000))=4 500 psia
T=70+2202=145 °F
Z = 0,95

Calculate Ims:

Ims=TZpms=605(.95)4 500=0,1277

Calculate pms:

pms=pts+0,01875γgHIms+Its=0,01875(.6)(10 000)0,1277+0,1140
pms = 4 000+465 = 4 465 psia

This is not close enough to the estimated value of 4 500 psia, therefore set


= 4 465 and repeat.

At T = 145, p = 4 465, Z = 0,95:

Ims=605(.95)4 465=0,1287
pms=4 000+112,50,1287+0,1140=4 000+464

pms = 4 464 psia which is close enough to the previously estimated value of 4 465.

Estimate pws:

  • pww*=pms(1+2,5·105(5 000))=4 464(1,125)=5 022 psia;
  • T = 220 °F;
  • Z = 1,027.

Calculate Iws:

Iws=680(1,027)5 022=0,1391

Calculate pws:

pws=pms+0,01875γgHIms+Iws=4 464+112,50,1287+0,1391
pws = 4 464+420 = 4 884 psia

For the second trial, Z = 1,019:

Iws=680(1,019)4 884=0,1419
pws=4 464+112,50,1287+0,1419=4 880 psia

A third iteration yields a value of pws = 4 880 psia. This compares with a value of 4 868 psia obtained using the average pressure and temperature method.

Flowing Bottom-Hole Pressure

For a flowing well the velocity is not zero, and ignoring acceleration, Equation 18 becomes, for a well inclined at an angle Ø from the vertical:

dpdL=ggcρ cos Ø +fρv22gcd.             Eq. 25

Several methods have been presented for integrating Equation 25 depending on the assumptions made for handling temperature and Z-factor. Only the average pressure and temperature and Cullender and Smith methods will be discussed.

Average Pressure and Temperature Method. Substituting the expression for gas density in terms of p, T, and Z into Equation 25 results in:

dpdL=pMZRTcos Ø +fv22gcd.           Eq. 26

Integration of Equation 26 assuming an average temperature in the flow string and evaluating Z at average conditions of pressure and temperature gives:

pwf2=ptf2 EXP(S)+25γgq2T¯Z¯f(MD)(EXP(S)1)Sd5,          Eq. 27


  • p = psia;
  • S


  • MD = measured depth, ft;
  • TVD = true vertical depth, ft;
  • T¯

    = °R;

  • q = MMscfd;
  • d = inches, and;
  • f = f(NRe, ∈/d) Jain or Colebrook equation or Picture 4.

The solution procedure is the same as for a shut-in well except for evaluation of the friction factor, which requires calculating a Reynolds number and estimating pipe roughness. Iteration is required since Z and μ must be evaluated at:


Dividing the well into several length increments and using the procedure described earlier will give more accurate results. Actually, any of the methods will give identical results if the well is divided into short enough increments.

Convergence is sometimes obtained faster if iteration is performed on the Z-factor rather than the unknown pressure. The procedure for this method is:

  1. Estimate Z* (A good first estimate is 0,9). Viscosity may be either estimated or evaluated at the known pressure.
  2. Calculate the unknown pressure using Equation 27 with Z = Z*.
  3. Calculate the average pressure,

  4. Evaluate Z and μ at



  5. Compare Z and Z*. If not close enough, set Z* = Z and go to Step 2. Repeat until abs(Z-Z*)/Z<0,001 or any other tolerance preferred. When the tolerance is met, the pressure calculated in Step 2 is the correct value.

Example 3 (a):

Use the average pressure and temperature method to calculate the flowing bottom-hole pressure for the following directional well:

  • γg = 0,75;
  • MD = 10 000 ft;
  • TVD = 7 000 ft;
  • Ts = 110 °F;
  • Tf = 245 °F;
  • ptf = 2 000 psia;
  • qsc = 4,915 MMscfd;
  • d = 2,441 in.;
  • = 0,0006 in.;
  • μ¯

    = 0,012 cp at p = 2 000 psia;

  • T




In terms of mass flow rate, the Reynolds number is:

NRθ=Cγgqscμd            Eq. 28


qsc = gas flow rateMMscfdMM m3/day
γg = gas gravity
μ = gas viscositycpkg/m-sec
d = pipe inside diameterin.m
C = constant20 01117,96
NRθ=20 011 γgqscμ¯d=20 011 (0,75)(4,915)0,012 (2,441)=2,518·106

From Equation 17 or Pic. 4, f = 0,015:

  1. Estimate Z* = 0,9:
    S=0,0375(0,75)(7 000)638 Z*=0,3086Z*.
  2. pwf2=(2 000)2 EXP (0,3086/Z*)+25(.75)(4,915)2(638)Z*(.015)(10 000)[EXP (0,3086/Z*)1]0,3086Z*(2,441)5
    pwf2=4·106 EXP (0,3086/Z*)+1,621·106(Z*)2 (EXP(.3086/Z*)1)
    For Z* = 0,9, pwf2

    = 5,636·106+536,966


    = 6,173·106, pwf = 2 485 psia.

  3. p¯=(ptf+pwf)/2=2 000+2 4852=2 242 psia.
  4. At p = 2 242 psia and T = 178 °F, Z = 0,806, μ = 0,012.
  5. abs(ZZ*)Z=0,90,8060,806=0,117, 

    which is too large.


For Z* = 0,806,


= 5,866·106+491,187


= 6,357·106, pwf

= 2 521 psia.


p¯=2 000+2 5212=2 261 psia.

4 At p = 2 261 psia and T = 178 °F, Z = 0,805.



which is close enough. Therefore pwf = 2 521 psia.

Cullender and Smith Method. Derivation of the Cullender and Smith method for flowing wells begins with Equation 26. The following substitutions are made for velocity:

v=qA,        q=qscpscTZTscpZsc,

which gives:

dpdL=pM cos ØZRT+MTZ psc2 f qsc2Rp Tsc2 2 gcd A2


pZTdpdL=MRpZT2 cos Ø+C,


C=8 psc2 qsc2 fTsc2 gc π2 d5

which is constant for a given flow rate in a particular pipe size. Separating the variables gives:

ptfpwfpZTdppZT2 cos Ø+C=MRoMDdL,        Eq. 29

which is applicable for any consistent set of units. Substituting field units and integrating the right-hand side of Equation 29 gives:

ptfpwfpZTdp0,001pZT2 TVDMD+F2=18,75γgMD,        Eq. 30


F2=0,667 f qsc2d5,


TVDMD=cos Ø.          Eq. 31

Writing Equation 30 in short notation and dividing the well into two increments of length H/2 gives:

Upper half of well:

18,75γg(MD) = (pmf-ptf)(Imf+Itf),

Lower half of well:

18,75γg(MD) = (pwf-pmf)(Iwf+Imf),


I=pTZ0,001pTZ2 TVDMD+F2             Eq. 32

The solution procedure is similar to that for the static case, but is more involved because of the more complicated definition of I. For practical purposes, F can be considered a constant since the only variable in the Reynolds number used in evaluating f is gas viscosity. Viscosity is a function of pressure, but for simplification of the calculations it can be evaluated at


and the known pressure.

Example 4:

The following data pertain to a flowing gas well. Use the Cullender and Smith method to calculate the flowing bottom-hole pressure.

  • γg = 0,75;
  • H = 10 000 ft;
  • Ts = 110 °F;
  • Tf = 245 °F;
  • ptf = 2 000 psia;
  • qsc = 4,915 MMscfd;
  • d = 2,441 in.;
  • = 0,0006 in.;
  • Ø = 0°;
  • μ¯

    = 0,012 cp (assume constant).


Calculate f and F2:

NRθ=20 011(0,75)(4,915)0,012(2,441)=2,518·106

From Equation 17: f = 0,015


Calculate Itf:

At p = 2 000 psia, T = 110 °F, Z = 0,71

pTZ=2 000(570)(.71)=4,942



(First Trial)

pmf*=2 000(1+2,5·105(5 000))=2 250 psia

Calculate Imf:

At p = 2 250, T = 110+67,5 = 178, Z = 0,797

pTZ=2 250(638)(.797)=4,425

Calculate pmf:

pmf=ptf+18,75γgHImf+Itf=2 000+18,75(.75)(10 000)197,81+181,60
pmf=2 000+371=2 371 (not close enough to pmf*)

Calculate Imf: (Second Trial)

At p = 2 371, T = 178, Z = 0,796

pTZ=2 371638(.796)=4,669

Calculate pmf:

pmf=2 000+140 625189,88+181,60=2 379

Calculate Imf: (Third Trial)

At p = 2 379, T = 178, Z = 0,796

pTZ=2 379638(.796)=4,684

Calculate pmf:

pmf=2 000+140 625189,41+181,60=2 379 psia

Therefore, the pressure at the mid-point of the well is 2 379 psia. The value of pwf is now calculated.



pwf*=2 379(1+2,5·105(5 000))=2 676

Calculate Iwf:

At p = 2 676, T = 245, Z = 0,867

pTZ=2 676705(0,867)=4,378

Calculate pwf: (First Trial)

pwf=pmf+140 625199,39+189,41=2 379+362=2 741

Calculate Iwf: (Second Trial)

At p = 2 741, T = 245, Z = 0,868

pTZ=2 741705(0,868)=4,479

Calculate pwf:

pwf=2 379+140 625196,00+189,41=2 744 psia

This is close enough to the previously calculated value of 2 741 psia. Therefore, the flowing bottom-hole pressure is 2744 psia.

Annular Flow. Many gas wells are dually completed, and one zone may produce through the tubing-casing annulus. This presents no problem in calculating static pressures, and either of the previously described methods can be used. For flowing wells either the Average Pressure and Temperature or the Cullender and Smith method may be used if the Hydraulic Radius concept is employed. The only modifications necessary are in calculating the effective diameter and area and the Reynolds number. It can be shown that the correct effective diameter is:

dh = dc-dt      Eq. 33


  • dh = effective diameter;
  • dc = casing inside diameter;
  • dt = tubing outside diameter.

This diameter is used in calculating both relative roughness and Reynold’s number. The effective area through which the gas is flowing is:


Therefore, when using either Equation 27 or the Cullender and Smith method, the d5 is replaced by:


There are no published data on values of roughness for an annulus. It would seem that a somewhat rougher surface than that for a round pipe would be applicable since the flow must pass the tubing collars. If one measurement of pressure drop and flow rate can be made, roughness can be back-calculated, and this value can be used for other flow rates.

The Cullender and Smith and Equation 27 methods can also be used to estimate pressure drop occurring during gas injection by using a negative value for F2 in Equation 30 or making the friction loss term negative in Equation 27.

Flow in Pipelines

For most practical cases of gas flow in pipelines, the line can be considered horizontal and the hydrostatic or elevation component, as well as the acceleration component, can be dropped from the general equation, Equation 18. This is not true if the line is transporting liquids, as will be discussed later. The general equation is:

dpdX=fρ v22gcd=pM f v2ZRT 2gcd            Eq. 34

Many solutions to Equation 34 have been proposed over the years. The difference in most of the solutions results from the method used to handle the friction factor, f, and the gas compressibility factor, Z. In most cases temperature is assumed constant and Z is evaluated at average pressure in the line. This requires an iterative solution if one of the pressures is unknown. Integration of Equation 34 over some distance, L, between upstream pressure p1 and downstream pressure p2 results in:

p12p22=25 γgq2T¯ Z¯ f Ld5,           Eq. 35


  • p = psia;
  • L = ft;
  • T¯

    = °R;

  • q = MMscfd (14,7 psia, 60 °F);
  • d = inches, and;
  • f = f(NRe, ∈/d) (Moody diagram or Jain eq.).

Equation 35 was derived using base or standard conditions as 14,7 psia and 60 °F. It can be put in a more general form by leaving the base conditions in the equation as variables. It is frequently advantageous to express Equation 35 in terms of flow rate; that is:

q=CTbpbp12p22γgf T¯ Z¯ L0,5 d2,5.            Eq. 36

The value of C depends on the units used in the equation. Table 1 gives the value of C for various combinations of units.

Example 5:

A pipeline is to deliver 320 MMscfd of gas to a downstream pressure of 600 psia. Using the following data, calculate the required upstream pressure.

  • T¯

    = 45 °F = 505 °R;

  • γg = 0,67;
  • d = 25,375 in.;
  • L = 100 miles = 528 000 ft;
  • = 0,0006 in.


Use Equation 35. p1 must be estimated to determine




. Estimate Δp/ΔL = 0,0005 psi/ft. Then:

p1*=p2+0,0005(528 000)=600+264=864 psia.
Table 1. Value of C for Various Units
psiaRin.ftscfd5 634
p¯=p*1+p22=864+6002=732 psia

At p = 732 psia, T = 45 °F, Z = 0,844, μ = 0,012 cp (assume constant).

Determine f:

NRθ=20 011 γgqμd=20 011(0,67)(320)0,012(25,375)=14,1·106

From Equation 17 (Jain Eq.), f = 0,0097

Calculate p1:

p12=(600)2+25(0,67) (320)2(505) (Z) (0,0097) (528 000) (25,375)5
p12=360 000+421 000 Z=360 000+355,898
p12=715,898        p1=846 psia

This is not close enough to the estimated value of 864. Set p*1 = 846 and re-calculate p1:

Second Trial:

p¯=600+8462=723 psia
Z = 0,846
p12=360 000+421 680(.846)=716 567
p1 = 847 psia

Since two successive values of p1, are approximately equal, the iterative solution has converged. The required upstream pressure is 847 psia.

Example 6:

Calculate the flow capacity of the following pipeline.

  • p1 = 4 140 kPa;
  • p2 = 2 760 kPa;
  • L = 30 000 m;
  • d = 0,152 m;
  • = 4,572·10-5m;
  • γg = 0,75;
  • Tb = 288 °K;
  • pb = 101,4 kPa;
  • T¯

    = 294 °K;

  • μ = 1,2·10-5kg/m-sec;
  • Z¯

    = 0,889.


The soiution will be trial-and-error or iterative since friction factor is required in Equation 36. This requires a value of Reynolds number, which requires flow rate, the unknown. Rather than estimate flow rate, it Is easier to estimate friction factor. The procedure is as follows:

  1. Estimate f*;
  2. Calculate q using Equation 36;
  3. Calculate N = f(q);
  4. Calculate f using the Moody diagram or Equation 17;
  5. Compare f and f*. If not close, set f* = f and go to Step 2. Continue until f = f*.
  1. Estimate f* = 0,01, a reasonable value for pipelines.
  2. q=1,149·106(288)101,4(4 140)2(2 760)20,75(f)(294)(0,889)(30 000)0,5(0,152)2,5
    q = 31 783/f0,5 = 31 783/(0,01) 0,5 = 317 830 m3
  3. NRθ=17,96γgqμd=17,96(0,75)q(1,2·105)(0,152)=7,385·106q
    N = 7,385·106(0,31783) = 2,347·106
  4. d=4,572·1050,152=0,0003
    f=1,142 log0,0003+21,25(2,347·106)0,92=0,015
  5. Comparing f and f*, 0,015 ≠ 0,010, therefore, set f* = 0,015.


q=31 783(0,015)0,5=259,507


N = 7,385·106(0,259507) = 1,916·106

4 f = 0,015

5 f = f*, therefore q = 259,507 m3/day.

The previous examples illustrate that if either the pressure drop or flow rate is unknown, the solution is iterative. If the unknown is diameter, the solution will also be iterative since diameter is necessary to evaluate the friction factor. This fact has prompted several investigators to substitute a specific equation for f into the general flow equation to make the solution for either q or d noniterative or explicit. The friction factor specified may be either diameter dependent only or Reynolds number dependent. None includes the dependence on pipe roughness.

The expressions incorporated into Equation 35 for f in several of the most popular pipeline equations are listed below.

Panhandle A0,085NRe0,147
Panhandle B0,015NRe0,183

Using these relationships for friction factor in Equation 36, the general form of the pipeline flow equation without f becomes:

q=a1ETbpba2p12p22T¯ Z¯ La3 1γga4da5         Eq. 37

where E is the efficiency factor, and the values of the ai constants used in the various equations are tabulated below.

The efficiency factor depends on the condition of the pipe and usually varies between about 0,7 and 0,92 for dry gas. It must be estimated just as pipe roughness was estimated in the previous equations.

Panhandle A435,871,07880,53940,46042,618
Panhandle B737,001,02000,51000,49002,530

The units to be used in Equation 37 are:

  • q = cu ft/day measured at Tb, pb;
  • T = R;
  • p = psia;
  • L = miles, and;
  • d = inches.

Example 7:

Using the following data, calculate the flow capacity of the pipeline using the Weymouth equation, the Panhandle B equation, and the Panhandle A equation.

  • p1 = 847 psia;
  • p2 = 600 psia;
  • d = 25,375 in.;
  • L = 100 miles;
  • γg = 0,67;
  • T¯

    = 505 °R;

  • Z¯

    = 0,846;

  • Tb = 520 °R;
  • pb = 14,7 psia;
  • E = 1,0.


p12p22T¯ Z¯ L=84726002505(0,846)(100)=8,366


q=433,5(35,374)1,0(8,366)0,510,670,5 (25,375)2,667
q = 301 610 000 cu ft/day = 301,6 MMscfd

Panhandle B:

q=737,0(35,374)1,02(8,366)0,51 10,670,49 (25,375)2,53
q = 359 732 857 cu ft/day = 359,7 MMscfd

Panhandle A:

q=435,87(35,374)1,0788 (8,366)0,5394 1.670,4604 (25,375)2,616
q = 364 247 375 cu ft/day = 364,2 MMscfd

The flow capacities using the various equations vary considerably and can be compared to the results of Example 5, which showed that the line would deliver 320 MMscfd at the specified pressures. If it is assumed that this is the correct value for q, efficiency factors of 1,061, 0,890, and 0,879 would have to be applied to the Weymouth, Panhandle B, and Panhandle A equations respectively to make the. calculated results agree.

The choice of which equation to use can be difficult. It is generally assumed that the Panhandle A equation is more applicable at Reynolds numbers in the transition region, and the B equation is more accurate in the fully turbulent region. If pipe roughness data are available, the value of friction factor can be determined and Equation 36 should be used.

All of the pipeline equations considered to this point have been applicable to horizontal or near-horizontal pipes only. That is, all of the calculated pressure drop is due to friction. If the gas changes elevation, an elevation component will apply even though it is very small for most pipelines. However, the elevation or hydrostatic component may be included in the general flow equation as follows:

q=CTbpb(p12p22 es) Sγgf T¯ Z¯ (es1)L0,5 d2,5,         Eq. 38


  • e = natural log base = 2,718;
  • s =

    0,0375 γgHZ¯ T¯,


  • H = elevation change, ft.

This method of handling the elevation change assumes that the rate of elevation change with pipeline length is constant. However, if this is not the case, very little error will result because of this assumption unless the average pressure in the line is very high, resulting in a large gas density effect.

Pipelines in Series

If a pipeline or gathering line consists of sections of different diameter pipe, the flow capacity of the entire pipeline can be calculated by first determining an equivalent length of a line of some arbitrary diameter that would have the same flow capacity as the system. The general flow equation or any of the specialized equations can then be used to determine the flow capacity. By using the fact that for a series pipeline:

pT=p1+p2+p3+=i=1N pi,

it can be shown that:

fLed5=i=1N fiLidi5,             Eq. 39

where N = number of pipe segments. If the specialized equations are used, the friction factor is not required and the equation becomes:

(Le)2a3(d)2a3=i=1N(Li)2a3(di)2a3,                  Eq. 40

where a3 and a5 are obtained from the previous table. For example, if the Weymouth equation is used, since a3 = 0,5 and a5 = 2,667, Equation 40 becomes:

Le=d5,333 i=1N Lidi5,333.              Eq. 41

Pipelines in Parallel

For pipes in parallel the total flow rate is the sum of the rates in the individual pipes, or:

qT=q1+q2+...+qN=i=1N qi.               Eq. 42

If the lengths of the individual pipes are the same, the total flow capacity is calculated from:

qT=Ci=1N di2,5fi0,5.         Eq. 43

When applied to the specialized equations, Equation 43 becomes:

qT=Ci=1N dia5.             Eq. 44

The diameter do of a single line having the same flow capacity as N parallel lines of diameters di is:

do=i=1N dia51/a5.             Eq. 45

In many cases only part of an existing pipeline will be paralleled or “looped” in order to increase its flow capacity. Equations 41, 44, and 45 can be combined to give:

qnew=qold1+Y1(1+W)210,5,           Eq. 46


  • qnew = new flow capacity after looping;
  • qold = flow capacity before looping;
  • Y = fraction of the original line that was paralleled starting at the outlet;
  • W =

    d2d12,5  f2f10,5 ;

  • d1 = diameter of original line;
  • d2 = diameter of parallel line;
  • f1 = friction factor for original line, and;
  • f2 = friction factor for parallel line.

Using the specialized equations, the expression for W becomes:

W=d2d1a5.             Eq. 47

Equation 46 can easily be solved for either W or Y to determine either the necessary line size for the parallel line or the length of line to be looped.

Effects of Liquids

The equations presented in the previous sections were derived for calculating the relationship between flow rate and pressure drop for dry gases. There are many cases in gas production operations in which some liquid will be traveling in the pipe along with the gas. These include gas wells producing some condensate or water and pipelines in which condensationmay occur or water may form. The presence of these liquids greatly increases the pressure drop for a given gas rate and reduces the efficiency of gathering systems.

If the liquid loading is low, the increased pressure drop in the tubing can be handled by adjusting the gas gravity used in the vertical gas flow equations. This will not suffice for pipeline flow because the liquids will accumulate in the low sections in the line. It then becomes necessary to apply a two-phase pressure drop method to design the piping system. Methods for use in both wells and pipelines will be discussed in this section.

Well Performance

The two-phase flow problem in flowing wells can be handled by either fluid gravity adjustment or by applying a two-phase flow correlation. There are many such correlations available, but only the Hagedorn and Brown method will be discussed in detail for vertical flow. The Beggs and Brill method, which will be used for pipeline flow, can also be used for well performance, but it has been found to sometimes overpredict pressure drop in vertical flow. Other well-known vertical, two-phase flow correlations are those of Poettmann and Carpenter, Orkiszewski, Duns and Ros and Gray.

Gravity Adjustment. The gravity adjustment procedure consists of adjusting the gas gravity to a mixture gravity to account for the added density due to the liquid and then using one of the equations presented earlier to calculate pressure drop in a well. Either Equation 27 or the Cullender and Smith method can be used. The mixture gravity is given by:

γm=γg+4 591 γL/R1+1 123/R,              Eq. 48


  • γm = adjusted fluid gravity (air = 1);
  • γg = dry gas gravity;
  • γL = liquid specific gravity, and;
  • R = producing gas-liquid ration, scf/STB.

The gravity adjustment method may be used with confidence for wells producing at high gas-liquid ratios. If a well is producing at a gas-liquid ratio of less than approximately 10 000 scf/STB, the two-phase correlations should be used. A gas-liquid ratio of 10 000 scf/STB expressed in terms of liquid loading is 100 STB/MMscf.

This simple method of adjusting the gas density for the presence of liquids assumes that the fluids are flowing as a homogeneous mixture. It should not be used for any gas-liquid ratio unless the well is flowing at a high enough rate to keep it from loading with liquids. The loading phenomenon is discussed in the next section.

Hagedorn and Brown Method. The Hagedorn and Brown method, ignoring acceleration, requires solution of Equation 25 for each increment that the well is divided into.

dpdh=ggcρm cos Ø +fρfvm22gcd            Eq. 49

Empirical correlations are presented for determining the mixture density ρm and the friction factor, f. The parameters in Equation 49 are defined by:

  • ρm = ρLHLg(1-HL);
  • ρL = liquid density;
  • ρg = gas density;
  • HL = liquid holdup (fraction of pipe occupied by liquid);
  • Ø = angle of well segment from vertical;
  • vm = vsL+vsg;
  • vsL = superficial liquid velocity = qL/Ap;
  • vsg = superficial gas velocity = qg/Ap;
  • Ap = area of flow string = πd2/4;
  • d = flow string ID;
  • ρf = ρn2m;
  • ρn = ρLλ+ρg(1-λ), and;
  • λ = vsL/vm.

The friction factor is calculated using the Jain equation or found from the Moody diagram using the pipe relative roughness and the following Reynolds number:

NRem=ρnvmdμm,             Eq. 50


  • μm =

    μLHL μg(1HL);

  • μL = liquid viscosity, and;
  • μg = gas viscosity.

Determination of HL requires the use of three empirical correlations. These are presented in Pictures 7, 8, and 9.

Scheme: correlation
Pic. 7 Correlation for viscosity number coefficient C
Scheme: holdup-factor
Pic. 8 Holdup-factor correlation
Scheme: second correction, correlation
Pic. 9 Correlation for secondary correction factor

In order to determine HL from these figures, the following dimensionless numbers must be evaluated from known data:

NLv = vsLL/gσ).25,
Ngv = vsgL/gσ).25,
Nd = d(ρL/gσ).5,


Nd = μL(g/ρL3).25,

where σ = gas-liquid surface tension. These equations are valid for any consistent set of units. For field units the equations are:

NLv = 1,938 vsLL/σ).25,
Ngv = 1,938 vsgL/σ).25,
Nd = 120,872 d(ρL/σ).5,


NL = 0,15726 μL(1,0/ρL3).25,


  • vsL, vsg = ft/sec;
  • ρL = lbm/cu ft;
  • σ = dynes/cm;
  • d = ft, and;
  • μL = centipoise.

The procedure for finding HL is:

  1. Calculate NL.
  2. Find CNL from Picture 7.
  3. Calculate

    where pa = base pressure (14,7 psia).

  4. Find


    from Picture 8.

  5. Calculate


  6. Find ψ from Picture 9.
  7. Calculate HL = ψ(HL/ψ).

A constraint on liquid holdup is that HL≥λ.

Once HL is determined, NRe and thus f can be calculated. The pressure gradient can then be calculated. This is Step 5 in the procedure for calculating a pressure traverse, which was presented earlier. All of the fluid properties and velocities used in the above equations are evaluated at the average pressure and temperature in the tubing increment.

When gas and liquid flow simultaneously in a pipe, the velocity of the mixture must be great enough so that the drag forces acting upward on the liquid are great enough to offset the gravity forces acting downward. Otherwise, the liquid will fall back and accumulate in the pipe. This increase in liquid holdup will increase the density of the mixture and therefore, the hydrostatic or elevation pressure gradient. The resulting increase in bottornhole pressure will result in a decrease in rate coming from the reservoir and therefore, a further decrease in velocity of the fluid in the pipe. If steps are not taken to increase the velocity, the well may die. The velocity of the fluid in the pipe depends on both the pipe size and the volumetric flow rate. The velocity may become too low if either the tubing is too large or the flow rate is too small. The liquid loading phenomenon is illustrated qualitatively in Picture 9 (a).

The liquid loading phenomenon
Pic. 9 (a) Effect of tubing size and flow rate on liquid loading

Most of the two-phase flow correlations mentioned earlier can be used to construct plots as illustrated in Picture 9 (a) and thus determine the minimum gas rate to avoid loading for a particular tubing size. A simple equation for predicting this minimum rate will be presented in a subsequent section.

Example 8:

During the calculation of a pressure traverse in a gas well producing liquid, the following conditions were determined at the average pressure and temperature in the pipe increment:

  • p = 1 500 psia;
  • T = 180 °F;
  • d = 2,992 in.;
  • μg = 0,012 cp;
  • μL = 0,45 cp;
  • σ = 25 dynes/cm;
  • vsg = 30 ft/sec;
  • vsL = 5 ft/sec;
  • = 0,0006 ft;
  • ρL = 50 lbm/cu ft;
  • ρg = 8 lbm/cu ft;

Using the Hagedorn and Brown method, determine the pressure gradient.


Before finding HL and f, some preliminary calculations are made:

vm = vsL+vsg = 5+30 = 35 ft/sec
λ = 5/35 = 0,143
ρπ = 50(0,143)+8(1-0,143) = 14 lb/ft3
Ap = πd2/4 = 0,7854(0,249)2 = 0,0487 ft2
ρL/σ = 50/25 = 2
NLv = 1,938(5)(2).25 = 11.52
Ngv = 1,938(30)(2).25 = 69,14
Nd = 120,872(0,249)(2).5 = 42,56
NL = 0,15726(0,45)[1/(50)(25)3].25 = 0,0024

Determine HL:

  1. NL = 0,0024;
  2. From Picture 7, CNL = 0,002;
  3. XH=11,52(0,002)(1 500)0,142,56(69,14)0,575 (14,7)0,1=7,53·105;
  4. From Picture 8, HL = 0,29;
  5. Xψ=69,14(0,0024)0,38(42,56)2,14=2,28·103;
  6. From Picture 9, ψ = 1,0;
  7. HL = 1,0(0,29) = 0,29;
ρm = 50(0,29)+8(1-0,29) = 20,18 lb/ft3
ρf = (14)2/20,18 = 9,71 lb/ft3
μm = (0,45)0,29 (0,012)(1-.29) = 0,034 cp
NRθm=1 488(14)(35)(0,249)0,034=5,29·106

From Picture 4 or Equation 17, f = 0,025:

dpdh=38,72 lb/ft3=0,269 psi/ft

A FORTRAN computer subroutine for the Hagedorn and Brown method is included in the appendix.

Pipeline Performance

Two methods will be presented for handling pipelines in which both liquid and gas are flowing. The Flanigan method accounts for the added pressure drop caused by lifting the liquid up the hills in a hilly-terrain pipeline and ignores any pressure recovery in the downhill sections. The angle of the hill is of no consequence in this method.

The Beggs and Brill correlation can be used to account for the hydrostatic or elevation pressure drop and is applicable to downward two-phase flow such as might occur in offshore gathering lines.

Flanigan Method. Flanigan proposed using the Panhandle A equation to calculate the pressure drop due to friction, based on the gas flow rate. A correlation for efficiency factor as a function of superficial gas velocity and liquid loading is presented in Picture 10.

Correlation for efficiency factor
Pic. 10 Flanigan efficiency factor

The gas velocity and liquid to gas ratio are in ft/sec and bbls/MMscf, respectively. The additional pressure drop due to hills is calculated from:

pe1=ρLHFΣh144               Eq. 51


  • Δpe1 = pressure drop due to hills, psi;
  • ρL = liquid phase density, lb/ft3;
  • HF = holdup factor, and;
  • Σh = sum of the vertical heights of the hills.

The holdup factor is a function of gas superficial velocity and is calculated from:

HF=11+0,3264 vsg1,006.             Eq. 52

The superficial gas velocity used in Picture 10 and Equation 52 must be calculated at the average pressure and temperature in the line. This requires an iterative solution since either p1 or p2 is unknown and


. If either q or d is unknown, they must be estimated before vsg can be calculated.

Beggs and Brill Method. The Beggs and Brill method requires the determination of the flow pattern that would exist in the pipeline if the pipe were horizontal. Different equations are used to calculate liquid holdup for each flow pattern. The flow patterns defined are shown in Picture 11.

Scheme: flow patterns
Pic. 11 Horizontal flow patterns

Determination of the correct flow pattern requires calculating several dimensionless numbers, including a two-phase Froude number.

The following variables are used to determine which flow pattern would exist if the pipe were in a horizontal position. This flow pattern is a correlating parameter and gives no information about the actual flow pattern unless the pipe is horizontal.

NFR=vm2g d
L1 = 316 λL0,302
L2 = 0,0009252 λL-2,4684
L3 = 0,10 λL-1,4516
L4 = 0,5 λL-6,738

The horizontal flow pattern limits are:

Limitsλ<0,01 and NFR<L1
orλL≥0,01 and NFR<L2
LimitsλL≥0,01 and L2<NFR≤L3
Limits0,01≤λL<0,4 and L3<NFR≤L1
orλL≥0,4 and L3<NFR≤L4
LimitsλL<0,4 and NFR≥L1
orλL≥0,4 and NFR>L4

When the flow falls in the transition region, the liquid holdup must be calculated using both the segregated and intermittent equations and interpolated using the following weighting factors.

HL(transition) = A·HL (segregated) + B·HL (intermittent),




B = 1-A.

The same equations are used to calculate liquid holdup for all flow patterns. The coefficients and exponents used in the equations are different for each flow pattern.

The liquid holdup depends on flow pattern and is calculated from:

HL(Ø) = HL(0)ψ,     Eq. 53

where HL(0) is the holdup that would exist at the same flow and pressure conditions in a horizontal pipe. It is calculated from:

HL(0)=a λLbNFRc,             Eq. 54

where a, b and c are determined for each flow pattern from Table 2.

Table 2
Flow Patternabc

The value calculated for HL(0) is constrained by:


The factor for correcting the holdup for the effect of pipe inclination is given by:

ψ = 1+C[sin (1,8Ø)-0,333 sin3 (1,8Ø)],

where Ø is the actual angle of the pipe from horizontal, and:

С=(1λL) ln (αλLe NLVf NFRg),

where α, e, f, and g are determined for each flow condition from Table 3.

Table 3
Flow Patternαθfg
Segregated uphill0,011-3,7683,539-1,614
Intermittent uphill2,960,305-0,44730,0978
Distributed uphillNo correctionC = 0, ψ = 1, HL = f(φ)
All flow patterns downhill4,70-0,36920,1244-0,5056

The value for C must be positive, and if a negative value is calculated, C is set equal to zero.

Once HL(Ø) is determined, the two-phase density is calculated from:

ρs = ρLHLgHg,

where Hg = 1-HL.

The pressure gradient due to elevation change is then:

dpdZe1=ggcρs sin Ø.                Eq. 55

The pressure gradient due to friction is:

dpdZf=ftpρnvm22 gcd,                 Eq. 56


ρn = ρLλLgλg

The no-slip friction factor fn is determined from the Moody diagram or from Equation 17 using the following Reynolds number:



μn = μLλLgλg.

The ratio of the two-phase to no-slip friction factor is calculated from:

ftpfn=es,             Eq. 57


S=[ln (y)]/{0,0523+3,182 ln (y)0,8725 [ln (y)]2+0,01853 [ln (y)]4},            Eq. 58



The value of S becomes unbounded at a point in the interval 1<y<1,2; for y in this interval, the function S is calculated from:

S = ln (2,2y-1,2).

Although the acceleration pressure gradient is very small except for high velocity flow, it should be included for more accurate results.

dpdZacc=ρsvmvsggcpdpdZ.             Eq. 59

If an acceleration term is defined as:


the total pressure gradient can be calculated from:

dpdZ=dpdZe1+dpdZf1Ek              Eq. 60

Example 9:

Using the following data for a hilly terrain pipeline, calculate the outlet pressure using the Beggs and Brill method.

  • qo

    = 7 140 STB/D;

  • qg

    = 25,7 MMcf/D;

  • d = 12 in.;
  • p1(inlet) = 425 psia;
  • T¯

    = 90 °F, γg = 0,70;

  • γo = 0,83 = 40 °API.

Divide the pipeline into two sections. Section 1 rises 300 ft in one mile. Section 2 drops 300 ft in 3 000 ft.


Section 1


Estimate Δp and calculate



p*=30 psi, p¯=42530/2=410 psia

2 From fluid property correlations, at 410 psia and 90 °F:

  • Rs = 96 scf/STB;
  • Bo = 1,047;
  • μo = 2,4 cp;
  • σo = 19,6 dyno/cm;
  • z = .925;
  • μg = 0,0105 cp.

3 Calculate flow rates and densities:

ρo=350 γo+0,0764 Rs γg5,615 Bo=350(.83)+.0764(96)(.7)5,615(1,047)
ρo = 50,29 lbm/cu ft
ρg=2,7 pγgZT=2,7(410)(.7).925(550)=1,52 lbm/cu ft
qo=6,49·105qoBo=6,49·105(7 140)(1,047)=0,485 ft3/sec
qg=3,27·107 Z(qgqo Rs)TP
qg=3,27·107(.925)[25,7·1067 140(96)](550) 410
qg = 10,1 ft3/sec

4 Calculate the in-situ superficial velocities:

vsL = qL/A = .485/.785(1)2 = 0,617 ft/sec
vsg = qg/A = 10,1/.785 = 12,87 ft/sec
vm = vsL+vsg = .617+12,87 = 13,49 ft/sec

5 Determine the flow pattern:

NFR = vm2/gd = (13,49)2/(32,2)(1) = 5,65
NLv=1,938 vsLρLσL.25=1,51
L1 = 316 λL0,302=124
L2 = 0,0009252 λL-2,4684=1,86
L3 = 0,10 λL-1,4516=8,82

Since λL>.01 and L2<NFR≤L3.

Flow pattern is transition, therefore interpolation is required.

6 Calculate the liquid holdup:

  1. Segregated:
    HL(0)=.98 λL.4846NFR.0888=.189
    C = (1- λL) ln (.011 λL-3,768 NLv3,539 NFR-1,614) = 5,516
    Ø = arcsin (300/5 280) = 3,257°
    ψ = 1+C [sin (1,8 Ø)-sin3 (1,8 Ø)/3] = 1,56
    HL(Ø)(segregated) = HL(0) ψ = .189 (1,56) = .295
  2. Intermittent:
    HL(0)=.845 λL.5351NFR.0173=.157
    C = (1-λL) ln (2,96 λL.305 NLv-.4437 NFR.0976) = .1246
    ψ = 1+C(.1015) = 1,013
    HL(Ø) (intermittent) = .157(1,013) = .159
    B = 1-A = .545
    HL (transition) = A·HL (seg.)+B·HL (int.) = .455(.295)+.545(.159) = .221

7 Calculate the actual and no-slip densities:

ρs = ρLHLgHg = 50,29(.221)+1,52(.779) = 12,3 lbm/cu ft
ρs = ρLλLgλg = 50,29(.0457)+1,52(.9543) = 3,75 lbm/cu ft

8 Calculate the friction factor:

NRθn=1 488 ρnvmdμLλL+μgλg=1 488(3,75)(13,49)(1)2,4(.0457)+.0105(.9543)=6,29·105

From Equation 17:

fn = .0126
X = ln y = -.066
S = X/[-0,0523+3,182X-0,8725X2+0,01853X4] = 0,248
ftp = fn EXP (S) = .0126(1,28) = 0,016

8 Calculate the pressure gradient:

dpdL=ggc ρθ sin Ø + ftpρnvm22 gcd1pθvmvsggcp=12,3(.057)+.016(3,75)(13,49)264,4(1)32,2(410)(144)=0,70+0,170,999=0,871 psf/ft=0,0061 psi/ft

10 Calculate the pressure drop:

p=dpdLL=.0061(5 280)=32,2 psi

For hand calculations this is ciose enough to the estimated Δp* of 30 psi. If a computer program were being used, the calculated Δp would be taken as the new estimated Δp* and the procedure would be repeated until the error was smaller. The pressure at the end of section 1 is 425-32,2 = 392,8 psia. The same procedure must be followed to calculate the pressure at the end of section 2.

Section 2


p*=0, p¯=392,8 psia;


  • Rs = 90;
  • Bo = 1,046;
  • μo = 2,4 cp;
  • μg = 0,0105 cp;
  • σo = 19,6 dyno/cm;
  • z = .925;


ρo=350(.83)+.0764(90)(.7)5,615(1,046)=50,28 lbm/cu ft

ρg=2,7(392,8)(.7).925(550)=1,46 lbm/cu ft
qo=.485(1,046)1,047=.484 ft3/sec
qg=3,27·107(.925)392,8·[25,7·1067 140(90)]550392,8=10,61 ft3/sec

4 vsL = .484/.785 = 0,616 ft/sec

vsg = 10,61/.785 = 13,52 ft/sec
vm = .616+13,52 = 14,14 ft/sec

5 λL = .616/14,14 = 0,4356

NFR = 14,142/(32,2) = 6,21
L1 = 316(.04356).302 = 122,7
L2 = 2,11    L3 = 9,45

Since λL>.01 and L2<NFR<L3, flow pattern is transition.


  1. Segregated:
    C = (1-λL) ln (4,7 λL-.3692
    Ø = -arcsin (300/3 000) = -5,74°
    ψ = 1+C (sin(1,8 Ø)-sin3(1,8 Ø)/3) = 1-1,75(.177) = .69
    HL(Ø) = .183(.69) = .126
  2. Intermittent:
    C = 1,75,   ψ = .69
    HL(Ø) = .153(.69) = .106
    A=9,456,219,452,11=.44   B=.56

7 ρs = 50,28(.115)+1,46(1-.115) = 7,07 lbm/cu ft

ρn = 50,28(.04356)+1,46(1-.04356) = 3,59 lbm/cu ft


NRθn=1 488(3,59)(14,14)2,4(.04356)+.0105(.9564)=6,59·105

fn=.0125     y=.04356(.115)2=3,29
X = ln      y = 1,192       S = .319
ftp = .0125 EXP (.319) = .017


dpdL=7,07(.10)+.017(3,59)(14,14)264,4(1).999=.707+.189.999=.518 psf/ft=.0036 psi/ft

10 Δp = (-.0036)(3 000) = -10,8 psi

The estimated Δp* was zero. Iteration through two more trials gives a pressure drop of -8,6 psi. The pressure at the outlet end of the pipe is then 425-32,2+8,6 = 401,4 psia.

Gas Flow Through Restrictions

There are several locations in the gas production system where the gas must pass through relatively short restrictions. Examples of these restrictions are perforations, subsurface safety valves, and surface chokes. The flow may be either critical or subcritical. In critical flow the velocity of the gas through the restriction is equal to the velocity of sound in the gas. Since pressure disturbances travel at the velocity of sound, a disturbance in the pressure downstream of the restriction cannot affect the upstream pressure or the flow rate.

In sub-critical flow the flow rate depends on both the upstream and downstream pressures. Surface chokes are usually sized so that flow will be critical, whereas in subsurface safety valves the flow is subcritical. Flow through well perforations will also be subcritical. Equations for both critical and subcritical flow are given in this section.

Flow through well perforations was discussed in Manual for Engineers to Calculate Amount of Gas in ReservoirsEssential Calculations of the Gas Reservoire Performance.

A general equation for flow through restrictions can be obtained by combining the Bernoulli equation with an equation of state and assuming that there are no irreversible or friction losses taking place. An empirical discharge coefficient is included to account for the simplifying assumptions used in deriving the equation. The following equation may be used for both critical (sonic) or subcritical (subsonic) flow. Tables 4 and 5 give values for the constants in the equation for various systems of units.

Table 4. Coefficients and Units for Equation 61
SymbolEnglish SystemMetric SystemSI Metric System
p abspsiakgf/cm2kPa
T abs°R°K°K
qsc=Cn(p1)(d)2γg(T1)Z1kk1p2p12/kp2p1k+1/k,           Eq. 61


  • qsc – volumetric gas flow rate;
  • Cn – coefficient based on system of units, discharge coefficient and standard conditions;
  • dID of bore opening to gas flow;
  • γg – gas specific gravity (air = 1,0), dimensionless;
  • k – ratio of specific heats = Cp/Cv, dimensionless;
  • p1 – upstream pressure, absolute units;
  • p2 – downstream pressure, absolute units;
  • T1 – upstream temperature, absolute units;
  • Z1 – compressibility factor at p1 and T1, dimensionless;
  • Cs – coefficient based on system of units;
  • Cd – discharge coefficient (empirical), dimensionless;
  • Tsc – standard temperature base, absolute units;
  • psc – standard pressure base, absolute units, and;
  • Rpc – critical pressure ratio, dimensionless.

Values of k can be obtained from:



  • M = molecular weight, lbm/mole;
  • Cp = specific heat, BTU/lbm – °R.

The pressure ratio at which flow becomes critical depends on the k value for the flowing gas and is given by:

Table 5. Coefficient for Equation 61
System of UnitsCdpsc absTsc absCn
English0,86514,696 psi491,68 °R799,06
0,86514,696 psi519,68 °R844,57
Metric0,8651,0332 kg/cm2273,16 °K371,83
0,8651,0332 kg/cm2288,72 °K393,01
SI Metric0,865101,325 kPa273,16 °K3,7915
0,865101,325 kPa288,72 °K4,0075

In calculating the values for Cn given in Table 5, a discharge coefficient of 0,865 was used. The discharge coefficient actually depends on the Reynolds number, the ratio of the diameter of the pipe to the diameter of the restriction, and the geometry of the restriction.

Example 10:

Using the following data, find the flow rate through the choke for:

a) p2 = 2 837 kPa, and;


  • p2 = 1 420 kPa;
  • d = 10 mm;
  • γg = 0,69;
  • k = 1,25;
  • Cd = 0,865;
  • psc = 101,325 kPa
  • Tsc = 288,72 °K;
  • T1 = 333 °K;
  • p1 = 3 546 kPa;
  • Z1 = 0,93.


From Table 5, Cn = 4,0075:



p2p1=2 8373 546=0.80,

, therefore flow is subcritical.

Using Equation 61:

qsc=4,0075(3 546)(10)2[(0,69)(333)(.93)]5   (5[(.8)1,6(.8)1,6])0,5
qsc = 97 213 (0,391) = 38 000 m3/d


p2p1=1 4203 546=0,4

, therefore flow is critical.

qsc = 97 213 (5[(.555)1,6-(.555)1,8])0,5
qsc = 97 213 (0,465) = 45 235 m3/d

Equation 61 has been modified for particular types of wellhead chokes. An equation which is used for the types of chokes manufactured by the Thornhill-Craver company is given below. This equation applies for 6 in. long chokes with rounded entrances operating in critical flow. The equation is:

qsc=605,4 A p1Cd(T γg)0,5               Eq. 61 (a)


  • qsc = flow rate, Mscfd;
  • A = area of choke opening, in.2;
  • p1 = upstream pressure, psia;
  • Cd = discharge coefficient, usually = 0,82;
  • T = upstream temperature, °R, and;
  • γg = gas specific gravity.

Example 11:

Recalculate the flow rate in Example 4·10b using Equation 4·61a.

  • d = 10 mm = 0,394 in.;
  • p1 = 3 546 kPa = 514 psia;
  • T1 = 333 °K = 600 °R;
  • γg = 0,69.


A = .7854(0,394)2 = 0,122 in.2
qsc=605,4(0,122)(514)(0,82)[(600)(0,69)]0,5=1 530 Mscfd=43 337 m3/d

An equation for calculating the pressure drop across a subsurface safety valve operating in subcritical flow was presented by the API in 1974. For the English System of units given in Table 4, the equation is:

p1p2=2,7 γgp1Z1T1(1β4) 6,23·104 Z1T1 qscp1d2 CdY,           Eq. 62


  • β = d/dp;
  • dp = pipe diameter, inches;
  • Cd = discharge coefficient (API suggests using 0,9), and;
  • Y = expansion factor;
  • Y=1[0,41+0,35β4] p1p2k p1.

The solution of Equation 62 is iterative since Y is a function of Δp = p1-p2. The value for Y ranges from about 0,67 to 1,0. For quick estimates of Δp, a value of 0,85 can be used.

Example 12:

A subsurface safety valve having a bean diameter of 1,0 in. is installed in a gas well equipped with 3,5 in. tubing (2,992 in. ID). The well is flowing at a rate of 20 MMscfd. Calculate the pressure drop across the SSSV if the pressure upstream of the SSSV is 2 000 psia. The temperature is 180 °F. Assume Cd = 0,9, Y = 0,85. Gas gravity is 0,70.


β=ddp=1,02,992=0,334      Z1=0,84
p=2,7(0,7)(2 000)0,84(640)(1(.334)4)·6,23·104(.84) (640) (20 000) (2 000) (1,0)2(.9) (.85) 
Δp = 6,994(4,378)2 = 134 psi

Use of Pressure Traverse Curves

Several equations were presented earlier for calculating flowing bottom-hole pressures in gas wells. Solution of these equations is iterative and, unless a computer is available, can require a considerable amount of time to solve. Estimates in the field can be made by using previously prepared pressure traverse curves, which can be calculated using conditions relative to a specific field. These curves can be prepared using either Equation 27 or the Cullender and Smith method. If small length increments are used in the calculations, Equation 27 is sufficiently accurate. Samples of flowing gas pressure traverse curves prepared using Cullender and Smith (see Pics. 12 through 16) are included for illustration.

Diagram: flowing gas gradients
Pic. 12 Vertical flowing gas gradients
Flowing gas gradients
Pic. 13 Vertical flowing gas gradients
Sample of flowing gas pressure
Pic. 14 Vertical flowing gas gradients
Diagram: flowing gas pressure
Pic. 15 Vertical flowing gas gradients
Scheme: gas pressure in tubes
Pic. 16 Vertical flowing gas gradients

An example problem illustrates the application procedure for determining flowing bottom-hole pressure. It should be pointed out that when using the prepared traverse curves, one has no control over gas gravity, flowing temperature, gas viscosity, or pipe roughness. These curves were taken from Brown.

Example 13:

Using the traverse curves, estimate the flowing bottom-hole pressure for the well producing under the following conditions:

  • H = 10 000 ft;
  • ptf = 1 200 psig;
  • d = 2,441 in.;
  • qsc = 9,8 MMscfd.


Select Picture 13 as the curve that closely matches the given conditions. The following procedure will yield an estimate of pwf:

  1. Starting at the known pressure of 1 200 psig on the pressure axis, draw a vertical line to intersect the 10 MMscfd traverse line.
  2. Draw a horizontal line to the depth axis. This locates the equivalent depth of 2 000 ft, which represents the wellhead.
  3. Add the well depth to the equivalent depth found in Step 2 to get 10 000 + 2000 = 12 000 ft.
  4. Draw a horizontal line from 12 000 ft to intersect the 10 MMscfd line.
  5. From this intersection, draw a vertical line to the pressure axis and read pwf = 2 150 psig.

Liquid Removal from Gas Wells

One of the most common factors that can reduce the deliverability of a gas well is the increased flowing wellbore pressure caused by the accumulation of liquids in the well bore. Liquid accumulation can occur in wells that have never produced large quantities of liquid.

When the velocity of the fluid in the tubing is too low to lift the liquids to the surface the liquid accumulates and begins to build up in the bottom of the well. The increased hydrostatic head acting on the formation further reduces the flow rate and thus the velocity of the fluid. This process continues until the well dies or begins to flow intermittently. The source of these liquids is either condensation of hydrocarbons or water from the gas or water entering from the formation.

Several methods have been used to either keep the liquid from accumulating or to remove the liquid as it accumulates. Some of these methods are discussed in this section.

Minimum Flow Rate for Continuous Liquid Removal

A model for calculating the minimum gas velocity for removing liquid droplets from wells was presented by Turner, et al. in 1969. The model was based on the fact that a freely falling particle in a fluid will reach a terminal velocity when the drag forces are equal to the gravitational forces. The terminal velocity is a function of the size, shape, and density of the particle and the density and viscosity of the fluid in which it is falling. The interfacial tension between the two fluids also influences this velocity since the droplet will shatter if the surface tension of the liquid is too small.

By writing a force balance on a droplet suspended in gas, the following equation was obtained:

vt=1,3 σ1/4(ρLρg)1/4Cd1/4 ρg1/2            Eq. 63


  • vt = terminal velocity of the droplet, ft/sec;
  • σ = interfacial tension, lbf/ft;
  • ρL = liquid density, lbm/ft3;
  • ρg = gas density, lbm/ft3, and;
  • Cd = drag coefficient.

The equation was adjusted using experimental data to evaluate the drag coefficient, and when surface tension is expressed in dynes/em, Equation 63 becomes:

vt=20,4σ1/4(ρLρg)1/4ρg1/2            Eq. 64

Turner, et al. also presented simplified versions of Equation 64 for field use that are easier to use since only the minimum pressure in the flow string is required. Equations were presented for use when the liquid is either water or condensate, which necessitate assuming average values for the fluid properties. The assumptions used were:

Liquidσ, dynes/cmρ, lbm/ft3

It was. also assumed that the gas density could be put in terms of pressure if a gas gravity of 0,6 and a temperature of 120 °F were used. This resulted in separate equations for water and condensate.

vg(water)=5,62(670,0031p)1/4(0,0031p)1/2        Eq. 65
vg(condensate)=4,02(450,0031p)1/4(0,0031p)1/2        Eq. 66

where the velocities are in ft/sec and the pressures are in psia. The minimum volumetric flow rate for a particular velocity and pipe size may be calculated from:

qsc (min.)=3,06vg ApT Z,             Eq. 67
  • qsc (min.) = minimum flow rate for continuous liquid removal, MMscfd;
  • vg = gas velocity, ft/sec;
  • A = conduit area, ft2;
  • T = flowing temperature, °R;
  • Z = gas compressibility factor evaluated at T and the pressure used to calculate vg, and;
  • p = wellhead pressure, psia.

Equations 65 and 66 were derived using wellhead flowing pressures since this was the only pressure measured in the field tests. Actually, the minimum velocity in a gas well will occur at the point of highest pressure, that is, at the bottom of the well. Using the bottom-hole flowing pressure to determine qsc (min.) will introduce a safety factor in the design.

A nomograph of Equation 67 is presented in Picture 17.

Gas rate calculating
Pic. 17 Nomograph for calculating gas rate required to lift liquids through tubing of various sizes (after Turner, et al.)

The product qZ is obtained from the nomograph, and to obtain qsc (min.), the Z factor must be calculated at surface conditions.

Example 14:

A gas well is producing against a fixed wellhead pressure of 1 150 psia through 3-1/2 in. (2,992 ID) tubing. Wellhead temperature is 140 °F and gas gravity is 0,70. Calculate the minimum flow rate required to keep this well unloaded il it produces salt water along with the gas.


vg=5,62(670,0031 p)0,25(0,0031 p)0,5
vg=5,62(670,0031 (1 150))0,25(0,0031 (1 150))0,5=8,4 ft/sec

at p = 1 150, T = 140 °F, γg = 0,7, the Z-factor is 0,86.

Aπ4 d2=.7854 (2,992/12)2=0,0488 ft2
qsc (min.)=3,06 p vg AT Z=3,06(1 150)(8,4)(0,0488)(600)(0,86)
qsc (min.) = 2,8 MMscfd

Turner, et al., observed that Ibeir melbod was applicable for liquid loading up to 130 bbl/MMscf. If bolb water and condensate are present, the equation for water should be used.

Liquid Removal Methods

Various melbods for removing liquids from gas wells, sometimes called dewatering, have been used in Ibe past. Descriptions of the methods and Ibe degree of success of Ibe various melbods are described in papers by Hutlas and Granberry and by Libson and Henry. Libson and Henry found Ibat Ibe minimum surface velocity for keeping wells unloaded in the Intermediate Shelf Area of soulbwest Texas was 1 000 ft/min or about 17 ft/sec. Melbods for liquid removal include:

  • pumping units;
  • gas lift;
  • plunger lift;
  • intermittent flow ilb flow controllers;
  • small tubing installation and;
  • soap injection.

Beam Pumping Units. Pumping units may be used to pump Ibe liquids up Ibe tubing, allowing the gas to be produced through the annulus. An advantage of using pumping units is Ibat they do not depend on gas velocity for lift and can be used to deplete Ibe field to a very low pressure. It is desirable to have the tubing set as close as possible to the bottom perforations or even below the perforations. A liquid cushion above the pump helps prevent gas from entering the pump.

On low liquid-rate wells, difficulties arise because a beam pumping unit cannot be adjusted for very low rates. The wells are sometimes pumped intermittently.

Plunger Lift. A plunger installed in a gas well producing liquid acts as an interface between the gas stored in the annulus and the liquid accumulated in the tubing above the plunger. A schematic of a plunger lift well is shown in Picture 18. The well remains shut in for a period of time and is then opened to allow the plunger, pushed by the annulus gas, to unload the liquids.

Scheme of a plunger lift
Pic. 18 Plunger lift operations

The automatic opening of the well may be accomplished by a motor valve on the flowline. The motor valve may be operated by a clock or by monitoring the flow rate of the well. If a clock is used the time cycles are adjusted by trial-and-error to find the optimum.

Surface flow controllers may be used that permit the well to flow until the gas velocity drops to some critical value. The well is then shut in for a time period that still must be determined by trial. The advantage of the flow controller is that it allows the well to flow for the maximum length of time before shut-in.

Small Tubing. Smaller tubing may be installed in gas wells as the flow rate decreases in order to maintain the velocity above a critical value. In some cases a smaller string may be run inside the existing tubing string, such as running 1 in. inside 2-7/8 in. tubing. This method works best on low volume wells in which friction loss is not severe.

Gas-Lift. A recent development in removing liquid from gas wells is the combination of a liquid diverter and gas-lift system. This system lifts liquid up the tubing string and produces gas through the annulus.

A liquid diverter is a device that opens when a predetermined head of liquid is accumulated above it, allowing the liquid to enter the tubing from the annulus. The tubing is open to atmospheric pressure at the surface. As liquid is diverted to the tubing, it accumulates until the gas-lift valve opens to lift the liquid slug to the surface.

The gas-lift valve may be actuated by the liquid head in the tubing if it is a so-called fluid operated valve. Various other means have been proposed to control the opening and closing of the gas-lift valve. Some of these are described by Hutlas and Granberr.

Soap Injection. Injection of surfactants or foaming agents into the annulus with a chemical pump and timer has produced successful results in some wells. By the reduction of the surface tension, water can be unloaded continuously in a foamed state. Soap injection has not been as successful in wells producing condensate because of the difficulty involved in obtaining a surfactant that will foam condensate.

Erosional Velocity

When fluid flows through a pipe at high velocities, erosion of the pipe can occur. This is especially true for high capacity gas flow in which the in-situ velocity may exceed 60 to 70 ft/sec. Erosion is not as much of a problem in oil wells, although some high gas-liquid ratio wells may be subject to erosion.

The velocity at which erosion begins to occur cannot be determined exactly, and if solid particles such as sand are in the fluid, erosion may occur at relatively low velocities.

The velocity at which erosion may occur has been related to the density of the fluid by the following equation.

ve=Cρ0,5              Eq. 68


  • ve = erosional velocity, ft/sec;
  • ρ = fluid density, lbm/ft3, and;
  • C = a constant that ranges between 75 and 150.

A good value for C has been found to be about 100. If C is set equal to 100, and the gas equation of state is used to express density, Equation 68 becomes:

ve=10029 p γgZRT0,5,

The equation may be expressed in terms of gas flow rate at standard conditions by:

qe=1,86·105 ApZT γg0,5,            Eq. 69


  • qe = erosional flow rate, Mscfd;
  • A = area of the pipe, ft2;
  • p = lowest pressure in the pipe, psia;
  • T = temperature at point where p is determined, °R;
  • Z = gas compressibility factor at p, T, and;
  • γg = gas gravity.

Example 15:

A gas well is producing through 2-7/8 in tubing at a wellhead pressure of 800 psia. The wellhead temperature is 140 °F and gas gravity is 0,65. Determine the maximum rate at which this well can produce without exceeding the erosional velocity.


A=π4d2=0,7854(2,441/12)2=0,032 ft2

At p = 800, T = 140 °F the value of Z is 0,91.

qθ = 1,86·105(0,032)(800/(.91)(600)(.65))0,5
qθ = 8 936 Mscfd = 8,9 MMscfd

Predicting Flowing Temperatures

All of the correlations presented previously require a value of fluid temperature in order to calculate the required fluid property or pressure drop. The flowing temperature profile in a gas well or an oil well is usually assumed to be linear between the surface temperature and the bottom-hole temperature. A linear temperature profile is also frequently assumed for surface flowline calculations. The linear assumption for well flow will usually not introduce significant errors if a good value for surface flowing temperature can be obtained. The heat loss from a fluid in a pipe is a function of the mass flow rate in the pipe, and will therefore change with a change in producing rate.

An algorithm for coupling pressure and heat loss calculations requires an iterative solution because the overall heat transfer coefficient and the enthalpy change depend on pressure. If some average heat transfer coefficient can be determined, an approximate temperature profile can be calculated independent of the pressure loss calculations. This will of course be less accurate, but in many cases the amount of data available will not be sufficient to perform the more accurate calculation.

Flowing Temperatures in Wells

An equation for temperature in a well as a function of location, L, as derived by Ramey, can be written as:

TL = T1-GT[L-A(1-EXP(-L/A))],     Eq. 70


  • T1 = temperature at fluid entry (L = 0);
  • TL = temperature at location L;
  • GT = geothermal gradient;
  • A = relaxation distance wCp/πdU;
  • w = mass flow rate, = ρscqsc;
  • Cp = specific heat of the flowing fluid;
  • d = pipe diameter;
  • U = overall heat transfer coefficient, and;
  • L = distance from fluid entry.

When the equation is written in this form it assumes that the fluid and surroundings temperature are equal at the inlet to the pipe. This will be the case for flowing wells where T1 is the reservoir temperature. Also included is the assumption that the heat loss is independent of time. This assumption limits application of Equation 70 to wells that have been producing for a considerable length of time.

When multiphase flow is occurring in a well, the variables involved in evaluating the relaxation distance, A, are very difficult to determine, especially the overall heat transfer coefficient U. In view of this fact, Shiu and Beggs developed an empirical method to estimate A based on measured temperature profiles from 270 wells. Using the measured temperatures, TL at various locations, L, a value of A for each test was calculated from Equation 70. An equation to estimate A was then developed as a function of data that will usually be known.

The equation is:

A=C1wC2ρLC3dC4(API)C5 γgC6,            Eq. 71


  • A = relaxation distance, ft;
  • w = total mass flow rate, lbm/sec;
  • ρL = liquid density at standard cond., lbm/ft3;
  • d = pipe ID, in.;
  • API = oil gravity, °API;
  • γg = gas gravity (air = 1);
  • C1 = 0,0149;
  • C2 = 0,5253;
  • C3 = 2,9303;
  • C4 = -0,2904, and;
  • C5 = 0,2608, and;
  • C6 = 4,4146.

Equation 71 is applicable for flowing oil wells only, although a similar approach could be used for gas wells.

Flowing Temperatures in Pipelines

In order to calculate a temperature profile in a pipeline it is usually assumed that the temperature of the surroundings is constant. Modification of Equation 70 to account for this results in:

TL = TS+(T1-TS) EXP (-L/A),     Eq. 72

where Ts is the surroundings temperature, and the other variables are defined in Equation 70.

For flow of gases the Joule-Thomson effect may be included, but since this effect depends on pressure, an iterative solution is required. The more rigorous equation is:

TL = TS+μA(dP/dL)+[T1-TS-μA(dP/dL)] EXP (-L/A),     Eq. 73


  • μ = Joule-Thomson coefficient, and;
  • dP/dL = pressure gradient at L.

If measured values of temperature and flow rate are available for one flow condition, a value for the group of terms Cp/πU may be calculated. The value of this group can be assumed constant for a particular well or even for an entire field if necessary. The procedure is:

  1. Using measured values of T and L, calculate A;
  2. Using the measured flow rate qsc, calculate
  3. For other conditions of qsc and/or d, calculate

This procedure can be used for both wells and pipelines.


Did you find mistake? Highlight and press CTRL+Enter

Декабрь, 28, 2022 88 0
Add a comment

Text copied