Diffusivity equation

Quick refresher

The equation is built from three parts. The first is the mass balance: what entered an element of the reservoir, minus what left, is what accumulated. The second is Darcy's law: the velocity is proportional to the pressure gradient (the rate at which pressure changes with distance). The third is compressibility: how much fluid the reservoir gives up when the pressure falls.
For radial inflow in a uniform reservoir, the “curvature” of the pressure profile equals the rate of change of pressure divided by the hydraulic diffusivity $\eta = k/(\phi\mu c_t)$ . This is the only reservoir parameter in the equation.
Several assumptions make the equation linear: the permeability, viscosity and compressibility do not depend on pressure, and the product of compressibility and pressure change is small, $c\,\Delta p \ll 1$ . For oil above the bubble point and for water the error is a fraction of a per cent. For gas the condition is not met.
There is one equation but many solutions. A particular solution is fixed by the initial condition and two boundary conditions: at the well (a constant rate or a constant bottomhole pressure) and at the outer boundary (an infinite reservoir, a closed boundary, or a constant pressure).
In the dimensionless variables $r_D$ , $t_D$ , $p_D$ no reservoir parameter is left in the equation. A solution found once is converted to any reservoir.
Linearity gives the main practical property: solutions can be added (the principle of superposition), and the pressure drop is proportional to the rate.

Key relationships

\frac{1}{r}\frac{\partial}{\partial r}\!\left(r\,\frac{\partial p}{\partial r}\right) = \frac{1}{\eta}\,\frac{\partial p}{\partial t},\qquad \eta = \frac{k}{\phi\mu c_t}

diffusivity equation for radial flow (SI)

\mathbf{u} = -\dfrac{k}{\mu}\left(\nabla p - \rho g\,\nabla D\right)

Darcy's law in vector form; $D$ is depth

r_D = \dfrac{r}{r_w},\quad t_D = \dfrac{k\,t}{\phi\mu c_t r_w^2},\quad p_D = \dfrac{2\pi kh\,(p_i - p)}{qB\mu}

dimensionless variables (SI)

t_D = 3.553\cdot10^{-3}\,\dfrac{k\,t}{\phi\mu c_t r_w^2}

metric oilfield units: mD, h, mPa·s, 1/MPa, m

\Delta p \approx \Delta p_{\text{lin}}\left(1 + \dfrac{c\,\Delta p_{\text{lin}}}{2}\right)

correction to the linear equation for a liquid of constant compressibility

Common pitfalls

Applying the equation to gas, and to oil below the bubble point, in the same way as to undersaturated oil. There the compressibility is large and changes with pressure.
Putting the compressibility of the oil alone into the diffusivity instead of the total $c_t$ , or the absolute permeability instead of the permeability to the mobile fluid.
Confusing the equation with its solution. A formula with the logarithm of time or with the term $\ln(r_e/r_w) - 3/4$ is a solution under certain boundary conditions, not a property of the equation itself.
Mixing units in the dimensionless time. The constant $3.553\cdot10^{-3}$ belongs to hours, mD, mPa·s, 1/MPa and metres, and $2.637\cdot10^{-4}$ belongs to hours, mD, cP, 1/psi and feet.

After this lesson you can

derive the diffusivity equation for radial flow from the mass balance, Darcy's law and the equation of state
list the assumptions, and explain what linearisation (making the equation linear) gives and when it is not allowed
write Darcy's law and the equation in general form, with a velocity vector, anisotropic permeability (permeability that depends on direction) and gravity
state the initial and boundary conditions of the main problems and link them to the flow regimes
convert to dimensionless variables (variables without units) and back, in metric and in oilfield units

Three laws, one equation

In the lesson Flow regimes the diffusivity equation appeared ready-made: on the left is the “curvature” of the pressure profile, on the right the rate at which it changes with time. Here we assemble it from its parts and see which assumption is made at each step. There are three parts.

Conservation of mass, written as the continuity equation. Fluid does not appear or vanish: the difference between the mass that entered an element of the reservoir and the mass that left it equals the change of mass inside the element.
Darcy's law, the equation of motion. The filtration velocity is proportional to the pressure gradient, and the coefficient of proportionality is the mobility $k/\mu$ (Rock permeability).
The equation of state. The density of the fluid and the pore volume depend on pressure, and compressibility sets this dependence (total compressibility).

The first statement is a general law of nature and always holds. The second and the third are models, and they carry the assumptions: slow (laminar, smooth) flow of one fluid, and a small, constant compressibility. So the limits of the equation are set by Darcy's law and the equation of state, not by the mass balance.

Analogy: a chain of capacitors and resistors

A reservoir can be pictured as a long electric chain. Each element of volume is a capacitor (it stores charge), and neighbouring elements are joined by resistors. Voltage corresponds to pressure, and current to flow rate. The resistance of a resistor corresponds to the viscous resistance of the rock, $\mu/k$ . The capacity of a capacitor corresponds to the storativity of the reservoir, $\phi c_t$ (its ability to store fluid by elastic expansion). If we start drawing current from one end of the chain, the voltage at first falls only on the nearest capacitors. To “reach” the far ones, the current must pass through all the resistors in between. The larger the resistance and the capacity, the slower the disturbance spreads; in a reservoir this means a small hydraulic diffusivity. In mathematics the analogy is exact: the voltage in such a chain obeys the same equation. In physics it is incomplete: in a reservoir it is matter that moves, and the properties of that matter themselves depend on pressure.

Darcy's law in general form

In the lesson on permeability, Darcy's law is written for flow along one axis. In a reservoir the fluid moves in three dimensions, and the filtration velocity is a vector: it has a size and a direction. Pressure, on the other hand, is a scalar field: one number at each point. What links them is the pressure gradient, a vector that shows in which direction the pressure grows fastest and how fast.

Darcy's law in vector form

\mathbf{u} = -\frac{k}{\mu}\left(\nabla p - \rho g\,\nabla D\right),\qquad \nabla p = \left(\frac{\partial p}{\partial x},\ \frac{\partial p}{\partial y},\ \frac{\partial p}{\partial z}\right)

where:

$\mathbf{u}$ filtration velocity vector: the flow rate through a unit of area perpendicular to the flow, m/s
$k$ permeability (here the same in all directions), m²
$\mu$ viscosity of the fluid, Pa·s
$\nabla p$ pressure gradient, Pa/m
$\rho$ density of the fluid, kg/m³
$g$ acceleration due to gravity, 9.81, m/s²
$D$ depth, counted downward; $\nabla D$ is a unit vector pointing down, m

Units: SI

Holds when

The flow of one fluid is slow (laminar). The irreducible water does not move, and $k$ is the permeability to the mobile fluid.
The rock is isotropic: the permeability is the same in all directions (the general case is below).
Inertial forces are neglected. For high-rate gas wells an inertial resistance is added to the viscous one (Inflow at high rates: non-Darcy flow).

The minus sign means that the fluid flows towards falling pressure. The term with $\rho g$ subtracts the hydrostatic part of the gradient. In a fluid at rest, the pressure grows with depth by $\rho g$ per metre, and this causes no motion. Only the departure from the hydrostatic distribution moves the fluid. If the $z$ axis points up, the bracket contains $\partial p/\partial z + \rho g$ . For flow along a bed that dips at an angle $\theta$ (the $x$ axis pointing down the dip), it contains $\partial p/\partial x - \rho g\sin\theta$ .

In calculations of inflow to a vertical well in a horizontal reservoir, the flow is horizontal and the gravity term drops out. It becomes important in dipping reservoirs, in crossflow between layers, and when oil is displaced by water or gas (Buckley–Leverett theory and the Welge method, Segregated flow and gravity displacement). There it is convenient to use the pressure brought to one reference depth, the datum depth, instead of the pressure itself (pressure brought to the datum depth).

The permeability of sedimentary rocks depends on direction (anisotropy). Suppose the coordinate axes point along the principal directions of permeability, usually along the bedding and across it. Then each component of the velocity is set by its “own” permeability: $u_x = -(k_x/\mu)\,\partial p/\partial x$ , and so on. But if the axes are turned relative to the principal directions (cross-bedding, or a model grid not aligned with the bedding), a pressure gradient along one axis also causes flow along the others. Then the permeability is written as a tensor, a table of nine numbers.

Darcy's law with a permeability tensor

\mathbf{u} = -\frac{1}{\mu}\,\mathbf{K}\,\nabla p,\qquad \mathbf{K} = \begin{pmatrix} k_{xx} & k_{xy} & k_{xz}\\ k_{yx} & k_{yy} & k_{yz}\\ k_{zx} & k_{zy} & k_{zz}\end{pmatrix}

where:

$\mathbf{K}$ permeability tensor, m²
$k_{xx},\ k_{yy},\ k_{zz}$ diagonal components: the link between flow along an axis and the gradient along the same axis, m²
$k_{xy},\ k_{xz},\ k_{yz}$ off-diagonal components: flow along one axis caused by the gradient along another, m²
$\mathbf{u},\ \mu,\ \nabla p$ filtration velocity vector, viscosity, pressure gradient, m/s, Pa·s, Pa/m

Units: SI

Holds when

Gravity is left out; it is added in the same way as in the isotropic case.
The tensor is symmetric ( $k_{xy} = k_{yx}$ and so on), so by turning the axes it can always be brought to a diagonal form with three principal permeabilities.

The product of the tensor and the vector is calculated like the product of a matrix and a column: $u_x$ contains $k_{xx}$ , $k_{xy}$ and $k_{xz}$ , multiplied by the three components of the gradient. The off-diagonal terms differ from zero only when the axes do not coincide with the principal directions. In analytical solutions for wells, a diagonal tensor or a single permeability is usually enough. A full tensor appears in reservoir simulation models (Grids in reservoir simulation models).

Mass balance in a ring of the reservoir

Take a uniform horizontal reservoir of thickness $h$ and a vertical well that penetrates its whole thickness. The flow is radial: all quantities depend only on the distance $r$ from the well axis and on time $t$ . Pick out in the reservoir a thin ring between the circles of radius $r$ and $r + \dd r$ .

SchematicA ring element of the reservoir (a sketch, not to scale): plan view and section. Fluid enters through the outer surface of the ring and leaves through the inner one towards the well. The difference goes into the change of the mass of fluid in the pores of the ring.

Through a cylindrical surface of radius $r$ , a mass flow rate $\rho u\cdot 2\pi rh$ passes towards the well, where $u$ is the filtration velocity directed towards the well. The amount that enters through the outer surface of the ring is not the same as the amount that leaves through the inner one. The difference equals the rate of change of the mass of fluid in the ring. Its pore volume is $2\pi rh\phi\,\dd r$ , and the mass of fluid in it is $\rho\phi\cdot 2\pi rh\,\dd r$ . Dividing the balance by $2\pi rh\,\dd r$ , we get the continuity equation:

Continuity equation and Darcy's law for radial flow

\frac{1}{r}\frac{\partial}{\partial r}\big(r\,\rho\,u\big) = \frac{\partial(\rho\phi)}{\partial t},\qquad u = \frac{k}{\mu}\,\frac{\partial p}{\partial r}

where:

$\rho$ density of the fluid at reservoir conditions, kg/m³
$u$ filtration velocity towards the well: the flow rate through a unit area of the cylindrical surface, m/s
$\phi$ porosity, fraction
$r,\ t$ distance from the well axis, and time, m, s
$k,\ \mu,\ p$ permeability, viscosity, pressure, m², Pa·s, Pa

Units: SI

Holds when

The flow is radial and horizontal, in a reservoir of constant thickness, and the well penetrates the whole thickness.
There is one mobile fluid, and Darcy's law applies.
The velocity $u$ is positive when the fluid flows towards the well. Pressure grows with distance, so there is a plus sign before $\partial p/\partial r$ .

The left side is the excess of the inflow into the ring over the outflow from it, per unit volume of the reservoir. The right side is the rate at which mass accumulates. In steady flow the right side is zero, the product $r\rho u$ does not depend on the radius, and the result is the Dupuit equation (Steady-state radial flow). When fluid is withdrawn from the reservoir, the right side is negative: the mass in the ring falls, and more leaves through the inner surface than enters through the outer one.

In three dimensions the same statement is written with the divergence, an operator that turns a flow vector into a number: how much mass “spreads out” of a unit volume in a unit time. The mass balance then takes the form $\nabla\cdot(\rho\mathbf{u}) = -\partial(\rho\phi)/\partial t$ . The sign is different here because $\mathbf{u}$ is a vector, not a velocity “towards the well”. The expression $\frac{1}{r}\frac{\partial}{\partial r}(r\,\ldots)$ is exactly the divergence of a radial flow.

Compressibility and linearisation

The continuity equation has two unknowns: pressure and density. The equation of state gives the link between them. For a liquid it is written with the compressibility $c = \frac{1}{\rho}\frac{\dd\rho}{\dd p}$ , the relative change of density per unit change of pressure. If the compressibility is constant, the density changes exponentially: $\rho = \rho_i\,\e^{c(p - p_i)}$ . The pore volume is described in the same way, by the pore compressibility $c_f = \frac{1}{\phi}\frac{\dd\phi}{\dd p}$ (Rock compressibility and the effect of stress on core properties). A fluid with a small, constant compressibility is called a slightly compressible fluid.

Accumulation of mass through the change of pressure

\frac{\partial(\rho\phi)}{\partial t} = \rho\phi\,(c + c_f)\,\frac{\partial p}{\partial t}\ \longrightarrow\ \rho\phi\,c_t\,\frac{\partial p}{\partial t}

where:

$c$ compressibility of the mobile fluid, 1/Pa
$c_f$ pore compressibility, 1/Pa
$c_t$ total compressibility $S_oc_o + S_wc_w + c_f$ , 1/Pa
$\rho,\ \phi,\ p,\ t$ density, porosity, pressure, time, kg/m³, —, Pa, s

Units: SI

Holds when

The compressibilities are constant over the pressure range considered.
The arrow marks the step to a reservoir whose pores hold, besides the mobile oil, immobile irreducible water. This water also expands when the pressure falls.

The right side of the continuity equation is now expressed through the pressure alone. How much mass a unit volume of the reservoir gives up when the pressure falls by one unit is set by the product $\phi c_t$ , the storativity of the reservoir (total compressibility). The permeability $k$ in such a reservoir is the permeability to oil at irreducible water saturation, not the absolute permeability.

Substituting Darcy's law and the equation of state into the continuity equation, and taking the permeability and viscosity as constant, we get an equation for the pressure alone, but a nonlinear one:

The equation before linearisation

\frac{\partial^2 p}{\partial r^2} + \frac{1}{r}\frac{\partial p}{\partial r} + c\left(\frac{\partial p}{\partial r}\right)^{2} = \frac{\phi\mu c_t}{k}\,\frac{\partial p}{\partial t}

where:

$p(r, t)$ pressure at the distance $r$ from the well axis at the time $t$ , Pa
$c$ compressibility of the mobile fluid, 1/Pa
$c_t$ total compressibility, 1/Pa
$\phi,\ \mu,\ k$ porosity, viscosity, permeability, —, Pa·s, m²

Units: SI

Holds when

The permeability and viscosity do not depend on pressure and are the same throughout the reservoir.
The compressibility of the fluid is constant.

The term with the square of the gradient appears because the density of the fluid changes along the flow. It makes the equation nonlinear: solutions cannot be added, and there are almost no ready-made formulas for it.

We compare the quadratic term with its neighbour $\frac{1}{r}\,\partial p/\partial r$ . Their ratio is $c\,r\,\partial p/\partial r$ . At the well, $r\,\partial p/\partial r$ is the “unit” of the pressure drop of radial flow, $qB\mu/(2\pi kh)$ . It is usually from a few tenths to several megapascals, while the compressibility of liquids is of the order of $10^{-3}$ MPa⁻¹. The ratio comes out at the order of thousandths, and the quadratic term is dropped. So the criterion is not that the gradient is small “in general” (the gradient has units, and at the well it is large). The criterion is that the product of the compressibility and the typical pressure change is small: $c\,\Delta p \ll 1$ .

Diffusivity equation for radial flow: the linear form and its assumptions

\frac{\partial^2 p}{\partial r^2} + \frac{1}{r}\frac{\partial p}{\partial r} = \frac{1}{r}\frac{\partial}{\partial r}\!\left(r\,\frac{\partial p}{\partial r}\right) = \frac{1}{\eta}\,\frac{\partial p}{\partial t},\qquad \eta = \frac{k}{\phi\mu c_t}

where:

$p(r, t)$ pressure at the distance $r$ from the well axis at the time $t$ , Pa
$\eta$ hydraulic diffusivity; in Russian-language literature it is also written $\kappa$ or $\chi$ , m²/s
$k$ permeability of the reservoir to the mobile fluid, m²
$\phi$ porosity, fraction
$\mu$ viscosity of the fluid, Pa·s
$c_t$ total compressibility, 1/Pa

Units: SI

Holds when

The reservoir is uniform and isotropic, of constant thickness. Porosity and permeability do not depend on pressure.
The flow is horizontal and radial: the well penetrates the whole thickness, and gravity is not taken into account.
There is one mobile fluid, and Darcy's law applies (flow without inertial losses).
The flow is isothermal (at constant temperature), and the viscosity is constant.
The compressibility is small and constant: $c\,\Delta p \ll 1$ .

The equation is linear: the pressure and its derivatives appear only to the first power, and the coefficient $\eta$ is constant. All the properties of the reservoir and the fluid are gathered in one parameter, the hydraulic diffusivity. The larger it is, the faster the pressure evens out: a high mobility $k/\mu$ speeds up the process, and a large storativity $\phi c_t$ slows it down.

In other unit systems (2)

Metric oilfield units (hydraulic diffusivity)

\eta\,[\text{m}^2/\text{h}] = 3.553\cdot10^{-3}\,\frac{k\,[\text{mD}]}{\phi\;\mu\,[\text{mPa·s}]\;c_t\,[\text{MPa}^{-1}]}

$k$ in mD, $\mu$ in mPa·s, $c_t$ in 1/MPa, $\eta$ in square metres per hour.

$3.553\cdot10^{-3} = 0.9869\cdot10^{-15}\cdot3600/(10^{-3}\cdot10^{-6})$ : it converts mD → m², hours → seconds, mPa·s → Pa·s and 1/MPa → 1/Pa.

Oilfield units

\eta\,[\text{ft}^2/\text{h}] = 2.637\cdot10^{-4}\,\frac{k\,[\text{mD}]}{\phi\;\mu\,[\text{cP}]\;c_t\,[\text{psi}^{-1}]}

$k$ in mD, $\mu$ in cP, $c_t$ in 1/psi, $\eta$ in square feet per hour.

$2.637\cdot10^{-4}$ comes from the same conversions, with pressure in psi (1 psi = 6894.76 Pa) and length in feet (1 ft = 0.3048 m).

CalculatedCalculated excess of the exact pressure drop over the drop given by the linear equation, against the product of the compressibility and the drop given by the linear equation,

c\,\Delta p

. The liquid has a constant compressibility, and the formula is in the advanced section. For oil and water

c\,\Delta p

rarely exceeds 0.02, and the error is no more than one per cent. For gas

c\,\Delta p

is of the order of

\Delta p/p

, that is, tenths; the curve gives only the order of magnitude for it.

Worked exampleLinearisation error for oil and for gas

Given (assumed data)

Well and reservoir (lesson Flow regimes): $q$ = 100 m³/d; $B$ = 1.3; $\mu$ = 1.2 mPa·s; $k$ = 50 mD; $h$ = 15 m
Oil compressibility $c$: $1.6\cdot10^{-3}$ MPa⁻¹
Fall of the bottomhole pressure after 64 h of work: 3.0 MPa
For comparison, gas: pressure 15 MPa, pressure drop 3.0 MPa

Solution

The “unit” of pressure drop: $\dfrac{qB\mu}{0.536\,kh} = \dfrac{156}{402} = 0.388$ MPa. The ratio of the dropped term to its neighbour at the well: $c\cdot0.388 = 0.00062$ .
Oil: $c\,\Delta p = 1.6\cdot10^{-3}\cdot3.0 = 0.0048$ . The exact solution for a liquid of constant compressibility gives a drop larger than the linear one by a factor of $-\ln(1 - x)/x = 1.0024$ .
That is 0.24%, or 0.007 MPa out of 3.0 MPa.
Gas: the compressibility is about $1/p = 0.067$ MPa⁻¹, so $c\,\Delta p \approx 0.2$ . The same formula gives a difference of about 12%.

Answer: For oil, the linear equation is wrong by 0.24% (0.007 MPa), much less than the error of a rate measurement. For gas, the difference is more than ten per cent.

What it means. So the formulas of the next lessons can be applied directly to oil and water. These are the formulas with the exponential integral and with the logarithm of time, and the Horner and Miller–Dyes–Hutchinson (MDH) methods. For gas we use them only after changing to pseudo-pressure (gas). The estimate for gas is only an order of magnitude: for gas, both the compressibility and the viscosity change with pressure. The data in this example are assumed for illustration.

When the linear equation applies

Fluid and rock: Water, oil above the bubble point
Product $c\,\Delta p$: 0.001–0.02
What is done: the linear equation; the error is up to 1%

Fluid and rock: Gas
Product $c\,\Delta p$: of the order of $\Delta p/p$ : 0.1–0.5
What is done: pseudo-pressure, and for large drops pseudo-time too

Fluid and rock: Oil below the bubble point
Product $c\,\Delta p$: large and changing
What is done: multiphase models (What a reservoir simulation model is and why it is needed)

Fluid and rock: Rock whose permeability falls with pressure
Product $c\,\Delta p$: —
What is done: pseudo-pressure with $k(p)$ , or a numerical model

Fluid and rock	Product $c\,\Delta p$	What is done
Water, oil above the bubble point	0.001–0.02	the linear equation; the error is up to 1%
Gas	of the order of $\Delta p/p$ : 0.1–0.5	pseudo-pressure, and for large drops pseudo-time too
Oil below the bubble point	large and changing	multiphase models (What a reservoir simulation model is and why it is needed)
Rock whose permeability falls with pressure	—	pseudo-pressure with $k(p)$ , or a numerical model

The values are approximate.

Initial and boundary conditions

The diffusivity equation is the same for all wells and reservoirs. What distinguishes the problems is what is given at the initial moment and at the boundaries of the area. For the solution to be unique, we need one initial condition and two boundary conditions: at the well and at the outer boundary.

Initial condition: Before the well starts, the pressure in the whole reservoir equals the initial pressure: $p(r, 0) = p_i$ . If the well has already been working, the initial state is the pressure distribution left by its history. We take it into account by superposition (Superposition principle).
Well: constant rate: The pressure gradient at the wall is given (formula below). This is the condition of well tests: the rate is held constant and the pressure is recorded. In English-language literature it is called constant terminal rate (constant terminal rate solution).
Well: constant bottomhole pressure: $p(r_w, t) = p_{wf}$ is given, and the rate falls with time (constant terminal pressure). This describes a well that works at a constant bottomhole pressure, and the inflow of water from an aquifer to an accumulation (Aquifer models: water influx).
Outer boundary: infinite reservoir: $p \to p_i$ as $r \to \infty$ . The condition suits any reservoir as long as the disturbance has not reached its boundaries.
Outer boundary: no flow: $\partial p/\partial r = 0$ at $r = r_e$ : a fault, a pinch-out, or the dividing line with a neighbouring well (no-flow boundary).
Outer boundary: constant pressure: $p(r_e, t) = p_i$ : fluid is fed in by injection or by an active aquifer.

The condition of a constant rate at the wall of the well

\left.r\,\frac{\partial p}{\partial r}\right|_{r = r_w} = \frac{qB\mu}{2\pi kh}

where:

$q$ rate at standard conditions; $qB$ is the rate at reservoir conditions, m³/s
$B$ formation volume factor of the fluid (at the initial reservoir pressure), m³/m³
$r_w$ radius of the well, m
$k,\ h,\ \mu$ permeability, thickness of the reservoir, viscosity, m², m, Pa·s

Units: SI

Holds when

The sandface rate (the rate at the bottom of the well) is constant from the start. The effect of the volume of the wellbore is not taken into account (Wellbore storage).
The inflow is uniform over the whole thickness of the reservoir.

This is Darcy's law written for the wall of the well: the whole rate passes through the surface $2\pi r_wh$ . The right side is the “unit” of the pressure drop of radial flow. It contains everything that makes one well different from another for the same dimensionless solution.

One equation, different solutions

Condition at the well: constant rate
Outer boundary: infinite reservoir
Solution: the exponential integral; at the well, the logarithm of time
Where it is covered: Line-source solution and the pressure function, Pressure drawdown analysis

Condition at the well: constant rate
Outer boundary: no flow
Solution: after the transient period, pseudo-steady-state flow
Where it is covered: Pseudo-steady-state flow and the shape factor

Condition at the well: constant rate
Outer boundary: constant pressure
Solution: in the limit, the Dupuit equation
Where it is covered: Steady-state radial flow

Condition at the well: constant pressure
Outer boundary: infinite or closed aquifer
Solution: falling rate and cumulative influx
Where it is covered: Unsteady-state water influx: the van Everdingen–Hurst method

Condition at the well: the rate changes; several wells
Outer boundary: any
Solution: a sum of the constant-rate solutions
Where it is covered: Superposition principle, Pressure build-up analysis

Condition at the well	Outer boundary	Solution	Where it is covered
constant rate	infinite reservoir	the exponential integral; at the well, the logarithm of time	Line-source solution and the pressure function, Pressure drawdown analysis
constant rate	no flow	after the transient period, pseudo-steady-state flow	Pseudo-steady-state flow and the shape factor
constant rate	constant pressure	in the limit, the Dupuit equation	Steady-state radial flow
constant pressure	infinite or closed aquifer	falling rate and cumulative influx	Unsteady-state water influx: the van Everdingen–Hurst method
the rate changes; several wells	any	a sum of the constant-rate solutions	Superposition principle, Pressure build-up analysis

CalculatedThree solutions of one equation for the well of the lesson Flow regimes (rate 100 m³/d, skin 0, circular area of 300 m radius). The initial condition and the condition at the well are the same, so the bottomhole pressure behaves in the same way until the disturbance reaches the outer boundary. Then it splits up depending on the condition there. The grey dashed line is the infinite reservoir. The time axis is logarithmic.

An important consequence follows: as long as the disturbance has not reached the boundaries, the solution for an infinite reservoir suits any real reservoir. So the same formula with the logarithm of time describes the start of every well's work, and the differences between reservoirs show up later.

The dimensionless form

In dimensional form the problem contains nine quantities: $k$ , $h$ , $\phi$ , $\mu$ , $c_t$ , $r_w$ , $q$ , $B$ and the initial pressure. There is no need to solve it again for each reservoir. It is enough to measure distances in well radii, pressure in “units” of $qB\mu/(2\pi kh)$ , and time in units of $r_w^2/\eta$ . This is the time in which a disturbance travels a distance of about one well radius. This gives the dimensionless variables.

Dimensionless variables: radius, time and pressure

r_D = \frac{r}{r_w},\qquad t_D = \frac{\eta\,t}{r_w^2} = \frac{k\,t}{\phi\mu c_t r_w^2},\qquad p_D = \frac{2\pi kh\,(p_i - p)}{qB\mu}

where:

$r_D,\ t_D,\ p_D$ dimensionless distance, time and pressure drop, —
$r_w$ radius of the well, m
$p_i$ initial reservoir pressure, Pa
$q,\ B$ rate at standard conditions and formation volume factor, m³/s, m³/m³
$k,\ h,\ \phi,\ \mu,\ c_t,\ \eta$ permeability, thickness, porosity, viscosity, total compressibility, hydraulic diffusivity, m², m, —, Pa·s, 1/Pa, m²/s

Units: SI

Holds when

The conditions of the diffusivity equation hold. The rate is constant, and the initial pressure is the same throughout the reservoir.
$p_D$ is positive when the pressure falls. For an injection well the sign is the opposite.

The dimensionless pressure shows how many “units” of the radial drop the pressure fall makes up. In steady-state inflow, at the wall of the well, it equals $\ln(r_e/r_w)$ (Steady-state radial flow). The dimensionless time compares the elapsed time with the characteristic time of the well $r_w^2/\eta$ . This is usually a fraction of a second or a few seconds, so in practical problems $t_D$ is measured in thousands and millions.

In other unit systems (2)

Metric oilfield units

t_D = 3.553\cdot10^{-3}\,\frac{k\,t}{\phi\mu c_t r_w^2},\qquad p_D = \frac{0.536\,kh\,(p_i - p)}{qB\mu}

$k$ in mD, $t$ in h, $\mu$ in mPa·s, $c_t$ in 1/MPa, $r_w$ and $h$ in m, $p$ in MPa, $q$ in m³/d at standard conditions.

$3.553\cdot10^{-3}$ is the same factor as for the hydraulic diffusivity; $0.536 = 2\pi\cdot0.0853$ is the constant of the Dupuit equation in these units.

Oilfield units

t_D = 2.637\cdot10^{-4}\,\frac{k\,t}{\phi\mu c_t r_w^2},\qquad p_D = \frac{kh\,(p_i - p)}{141.2\,qB\mu}

$k$ in mD, $t$ in h, $\mu$ in cP, $c_t$ in 1/psi, $r_w$ and $h$ in feet, $p$ in psi, $q$ in STB/d (barrels of oil at standard conditions per day).

$141.2 = 1/(2\pi\cdot1.127\cdot10^{-3})$ , where $1.127\cdot10^{-3}$ is the constant of Darcy's law in oilfield units.

The dimensionless problem of a well at a constant rate

\begin{aligned} &\frac{1}{r_D}\frac{\partial}{\partial r_D}\!\left(r_D\,\frac{\partial p_D}{\partial r_D}\right) = \frac{\partial p_D}{\partial t_D}\\[4pt] &p_D = 0\ \text{at}\ t_D = 0,\qquad r_D\,\frac{\partial p_D}{\partial r_D} = -1\ \text{at}\ r_D = 1\end{aligned}

where:

$r_D,\ t_D,\ p_D$ dimensionless distance, time and pressure drop, —

Units: dimensionless

Holds when

The same as for the dimensional equation; a constant rate.
The third condition is set at the outer boundary: $p_D \to 0$ as $r_D \to \infty$ , or $\partial p_D/\partial r_D = 0$ , or $p_D = 0$ at $r_D = r_{eD}$ , where $r_{eD} = r_e/r_w$ .

No reservoir parameter is left, either in the equation or in the condition at the well. The solution $p_D(r_D, t_D)$ for an infinite reservoir is one function of two variables, common to all reservoirs. For a bounded reservoir, one more parameter appears: the dimensionless radius $r_{eD}$ .

The way back is simple: $p = p_i - \dfrac{qB\mu}{2\pi kh}\,p_D$ . Other dimensionless times are also convenient: those related not to the well radius but to the size of the area, $t_D/r_{eD}^2 = \eta t/r_e^2$ and $t_{DA} = \eta t/A$ . They tell us whether the disturbance has reached the boundaries (time to reach the regime).

Worked exampleDimensionless variables of the well in metric and oilfield units

Given (assumed data)

Reservoir and oil (lesson Flow regimes): $k$ = 50 mD; $\phi$ = 0.20; $\mu$ = 1.2 mPa·s; $c_t = 2.1125\cdot10^{-3}$ MPa⁻¹
Well: $r_w$ = 0.1 m; $h$ = 15 m; $q$ = 100 m³/d; $B$ = 1.3
Unit conversions: 1 m = 3.2808 ft; 1 m³ = 6.2898 bbl; 1 MPa = 145.04 psi

Solution

Metric units. The denominator: $\phi\mu c_t r_w^2 = 5.07\cdot10^{-6}$ . After 1 h: $t_D = \dfrac{3.553\cdot10^{-3}\cdot50}{5.07\cdot10^{-6}} = 35{,}039$ .
The pressure “unit”: $\dfrac{qB\mu}{0.536\,kh} = \dfrac{156}{402} = 0.388$ MPa.
Oilfield units: $c_t = 1.4565\cdot10^{-5}$ psi⁻¹, $r_w = 0.3281$ ft, $h = 49.21$ ft, $q = 629$ STB/d.
The denominator: $\phi\mu c_t r_w^2 = 3.7625\cdot10^{-7}$ . After 1 h: $t_D = \dfrac{2.637\cdot10^{-4}\cdot50}{3.7625\cdot10^{-7}} = 35{,}043$ .
The pressure “unit”: $\dfrac{141.2\,qB\mu}{kh} = 56.3$ psi, which is 0.388 MPa.
After 6 h of work, $t_D = 210{,}237$ . At the wall of the well in an infinite reservoir, $p_D = \tfrac12(\ln t_D + 0.809) = 6.533$ (Line-source solution and the pressure function), so the pressure fall is $6.533\cdot0.388 = 2.535$ MPa.

Answer: In both systems of units

t_D \approx 35{,}039

per hour, and the pressure “unit” is 0.388 MPa (56.3 psi). The dimensionless quantities do not depend on the choice of units.

What it means. The agreement is a check of the constants: $3.553\cdot10^{-3}$ and $2.637\cdot10^{-4}$ , 0.536 and 141.2 describe the same thing in different units. A difference in the fourth or fifth significant digit is a trace of the rounding of the constants and of the conversion factors. The characteristic time of the well $r_w^2/\eta$ is $1/35{,}039$ h, which is about 0.1 s. After only a minute of work the dimensionless time is already measured in hundreds. The data in this example are assumed for illustration.

Worked exampleOne calculation, two reservoirs

Given (assumed data)

Reservoir A: the well of the previous example: $k$ = 50 mD, $\mu$ = 1.2 mPa·s
Reservoir B: $k$ = 10 mD, $\mu$ = 2.4 mPa·s; the other data are the same
Rate: 100 m³/d in both cases

Solution

The mobility $k/\mu$ in reservoir B is smaller by a factor of $5\cdot2 = 10$ . The hydraulic diffusivity is smaller by the same factor, so the same $t_D$ corresponds to a time 10 times longer.
What takes 1 h in reservoir A takes 10 h in reservoir B. The dimensionless time is 35,039 in both cases, and the dimensionless pressure at the well is $p_D = 5.637$ .
The pressure “unit” is proportional to $\mu/k$ : in reservoir A it is 0.388 MPa, in reservoir B 3.88 MPa.
The fall of the bottomhole pressure is 2.19 MPa in reservoir A after 1 h, and 21.9 MPa in reservoir B after 10 h.
For the pressure fall in reservoir B to stay the same as in reservoir A, the rate must be cut by a factor of 10, to 10 m³/d.

Answer: A tenfold fall in mobility stretches all times by a factor of 10 and increases the pressure drop at the same rate by the same factor.

What it means. The dimensionless solution is one, but each reservoir has its own scales: time is proportional to $\phi\mu c_tr_w^2/k$ , and pressure to $qB\mu/(kh)$ . A drop of 21.9 MPa is comparable to the reservoir pressure itself, which means reservoir B cannot deliver 100 m³/d. So the dimensionless solution shows at once what rate and what test duration can be planned in a new reservoir. The data in this example are assumed for illustration.

What the equation tells us about pressure before we solve it

A lot about the behaviour of the pressure can be seen from the equation itself, without the formulas of the solutions.

A disturbance spreads out; it does not run like a wave. The equation contains the first derivative with respect to time and the second with respect to distance. So the distance covered by the disturbance grows as $\sqrt{\eta t}$ : to “see” twice as far, we need four times the time (radius of investigation).
Solutions can be added. The equation is linear, so the sum of two solutions is also a solution. The pressure from several wells, or from a well whose rate changes, equals the sum of the pressures from simple components (Superposition principle). The analysis of pressure build-up is built on this (Pressure build-up analysis).
The drop is proportional to the rate. Doubling the rate doubles the pressure fall at every point and at every moment. So it is convenient to express test results through the ratio of the pressure drop to the rate.
Near the well the profile is almost steady-state. Where the disturbance passed long ago, the right side of the equation is small compared with each term on the left. So the pressure depends linearly on $\ln r$ , as in the Dupuit equation. The transient character shows at the edge of the drawdown cone (the region where the pressure has dropped).
Sharp differences are smoothed out. Any unevenness of the pressure profile spreads out with time, and the reverse never happens. So from a pressure record in a well we cannot recover unambiguously the small details of the reservoir structure far from it (Well testing: purposes and types of test).
The speed of spreading is formally infinite. According to the equation, the pressure at any point starts to change right after the well starts. But far away this change is vanishingly small, so we speak of a finite radius of investigation.

How exactly the pressure depends on distance and time at a constant rate is shown by the line-source solution, the subject of the next lesson (Line-source solution and the pressure function).

Terms and symbols

Continuity equation: A statement of conservation of mass: the difference between the inflow into an element of the reservoir and the outflow from it equals the rate at which mass accumulates in the element.
Equation of state: The link between the density of the fluid (and the pore volume) and the pressure. For a liquid it is set by the compressibility.
Hydraulic diffusivity $\eta$: $k/(\phi\mu c_t)$ , in m²/s. It is the only reservoir parameter in the equation and sets the speed at which the pressure evens out.
Linearisation: Dropping the term with the square of the pressure gradient and replacing the properties of the fluid and the rock by constants. It is allowed when $c\,\Delta p \ll 1$ .
Initial and boundary conditions: The pressure at the initial moment, the condition at the well (a rate or a pressure) and the condition at the outer boundary. Together with the equation, they fix a single solution.
Dimensionless variables $r_D$ , $t_D$ , $p_D$: Distance in well radii, time in units of $r_w^2/\eta$ , and pressure fall in units of $qB\mu/(2\pi kh)$ .
Gradient and divergence: The gradient turns the pressure field into a vector of its steepest growth. The divergence turns a flow vector into a number: how much fluid “spreads out” of a unit volume.

Common mistakes and misconceptions

Wrong

“The diffusivity equation is just Darcy's law written in another way.”

In fact

Darcy's law is one of three components. Without the mass balance and the compressibility, only steady-state flow follows from it. Time enters the problem thanks to the compressibility.

Wrong

“The square of the gradient is neglected because the pressure gradient is small.”

In fact

At the well the gradient is in fact large. What is small is the product of the compressibility and the pressure change, $c\,\Delta p$ : for oil and water it is thousandths and hundredths.

Wrong

“There is one equation, so there is one pressure formula for all wells.”

In fact

The solution is also decided by the initial and boundary conditions. An infinite reservoir, a closed area and an area fed through its boundary give three different formulas. Only their early part is common.

Wrong

“The oil compressibility goes into the hydraulic diffusivity.”

In fact

It is the total compressibility $c_t$ that goes in: that of the oil, the irreducible water and the pores. In the example of the lesson on flow regimes, the pore compressibility gave 38% of the total.

Wrong

“Dimensionless time is just time in other units.”

In fact

The time scale $r_w^2/\eta$ is different for each reservoir and each well. In a reservoir with ten times lower mobility, the same $t_D$ corresponds to a time ten times longer.

In depth: advanced detail

How accurate the linearisation is: the equation for density

We do not have to drop the quadratic term; we can remove it by changing the variable. For a liquid of constant compressibility $\rho\,\partial p/\partial r = \frac{1}{c}\,\partial\rho/\partial r$ . So the mass flux is proportional to the gradient of the density, and the continuity equation turns into an equation for the density, which is linear with no approximation at all:

Diffusivity equation for density

\frac{1}{r}\frac{\partial}{\partial r}\!\left(r\,\frac{\partial\rho}{\partial r}\right) = \frac{\phi\mu c_t}{k}\,\frac{\partial\rho}{\partial t}

where:

$\rho(r, t)$ density of the fluid at reservoir conditions, kg/m³
$\phi,\ \mu,\ c_t,\ k$ porosity, viscosity, total compressibility, permeability, —, Pa·s, 1/Pa, m²

Units: SI

Holds when

There is one mobile fluid with a constant compressibility $c$ . The permeability, viscosity and the product $\phi c_t$ are constant.
The withdrawal is given as a constant mass rate (a constant rate at standard conditions). Then the condition at the well is also linear in the density.

The density obeys the same equation as the pressure in the linear theory, with the same coefficient and the same conditions. So the formulas of the linear theory describe the relative change of the density exactly: $\rho/\rho_i = 1 - c\,\Delta p_{\text{lin}}$ , where $\Delta p_{\text{lin}}$ is the pressure drop from the linear equation.

The exact pressure is found from the density through the equation of state $\rho = \rho_i\,\e^{-c\,\Delta p}$ :

The exact pressure drop through the drop from the linear equation

\Delta p = -\frac{1}{c}\ln\left(1 - c\,\Delta p_{\text{lin}}\right) \approx \Delta p_{\text{lin}}\left(1 + \frac{c\,\Delta p_{\text{lin}}}{2}\right)

where:

$\Delta p$ exact pressure fall $p_i - p$ at a given point at a given moment, Pa
$\Delta p_{\text{lin}}$ the same from the linear diffusivity equation, Pa
$c$ compressibility of the mobile fluid, 1/Pa

Units: any consistent units (the product $c\,\Delta p$ has no dimension)

Holds when

The same as for the density equation; the approximate form holds for $c\,\Delta p_{\text{lin}} \ll 1$ .

The linear equation underestimates the pressure drop, and the relative error is about half of the product $c\,\Delta p$ : at 0.01 it is 0.5%, and at 0.1 about 5%. The plot in the main section is built from this formula.

General form: three dimensions and anisotropy

If we put Darcy's law in vector form into the mass balance and linearise in the same way as for radial flow, we get the diffusivity equation for any geometry:

Diffusivity equation in three dimensions

\nabla\cdot\left(\frac{\mathbf{K}}{\mu}\,\nabla p\right) = \phi c_t\,\frac{\partial p}{\partial t}\quad\Rightarrow\quad k_x\frac{\partial^2 p}{\partial x^2} + k_y\frac{\partial^2 p}{\partial y^2} + k_z\frac{\partial^2 p}{\partial z^2} = \phi\mu c_t\,\frac{\partial p}{\partial t}

where:

$\mathbf{K}$ permeability tensor, m²
$k_x,\ k_y,\ k_z$ principal permeabilities (the coordinate axes point along the principal directions), m²
$p(x, y, z, t)$ pressure, Pa
$\phi,\ \mu,\ c_t$ porosity, viscosity, total compressibility, —, Pa·s, 1/Pa

Units: SI

Holds when

The equation is linearised. Gravity is not taken into account (or the pressure is replaced by the pressure brought to the datum depth).
The second form is for a uniform reservoir with the axes along the principal directions of permeability.

For an isotropic reservoir, the left side becomes $k\,\nabla^2 p$ , where $\nabla^2$ is the Laplace operator, the sum of the second derivatives with respect to the three coordinates. In a reservoir simulator this equation is written for each grid cell and solved numerically (What a reservoir simulation model is and why it is needed).

Anisotropy in the plan view is removed by stretching the coordinates. If we put $x' = x\sqrt{\bar k/k_x}$ and $y' = y\sqrt{\bar k/k_y}$ , where $\bar k = \sqrt{k_xk_y}$ , the equation for planar flow takes an isotropic form with permeability $\bar k$ . So a test of a vertical well in a reservoir with different permeabilities in the two horizontal directions gives the geometric mean $\sqrt{k_xk_y}$ . In plan, the drawdown cone is not a circle but an ellipse, stretched along the direction of higher permeability.

When the equation stops working

Gas. The compressibility is of the order of $1/p$ , and the viscosity depends on pressure. The pressure is replaced by the real gas pseudo-pressure, after which the equation is almost linear again (gas in the flow regimes, pseudo-pressure).
Oil below the bubble point. Two phases move in the pores, the mobility of each depends on the saturation (Relative permeability), and the compressibility of the system rises sharply. We need equations for each phase, and they are solved numerically.
Permeability that depends on pressure. In fractured and weakly cemented rocks the permeability falls when the pore pressure falls (Rock compressibility and the effect of stress on core properties). The equation becomes nonlinear.
High velocities near the wellbore. Darcy's law is supplemented with an inertial term. In well solutions it shows up as a skin that grows with the rate (Inflow at high rates: non-Darcy flow).
Fractured-porous reservoirs. The fractures conduct flow, and the rock blocks store the fluid. The pressure in the fractures and in the blocks is described by two coupled diffusivity equations.

Connections to other topics

Flow regimes — the three flow regimes as three kinds of solutions of the equation; total compressibility.
Rock permeability — Darcy's law and the units of permeability.
Line-source solution and the pressure function — the solution for a well at a constant rate in an infinite reservoir.
Steady-state radial flow — the solution with a zero right side: the Dupuit equation.
Pseudo-steady-state flow and the shape factor — the solution for a closed area after the transient period.
Superposition principle — adding solutions, a consequence of the linearity of the equation.
Well testing: purposes and types of test — how the solutions of the equation are used to test wells.
Rock compressibility and the effect of stress on core properties — pore compressibility and the dependence of rock properties on stress.
What a reservoir simulation model is and why it is needed — the numerical solution of the same equation on a grid of cells.

Check yourself

1
Which of the statements on which the diffusivity equation is built brings dependence on time into the problem?
1. A.The constancy of the viscosity
2. B.The mass balance together with the compressibility of the fluid and the pores
3. C.The condition of a constant rate at the well
4. D.Darcy's law
Show the answer
Correct answer: B. The mass balance together with the compressibility of the fluid and the pores
Darcy's law links the velocity with the pressure gradient at a given moment and contains no time. Time appears on the right side of the mass balance, in the rate at which fluid accumulates. Fluid can accumulate or decrease only because the fluid and the pores are compressible. In an incompressible system the right side is zero, and steady-state flow is what remains.
2
A reservoir has twice the permeability of its neighbour, but three times the total compressibility (there is free gas in the pores). The porosity and viscosity are the same. How do the hydraulic diffusivities compare?
1. A.They are the same
2. B.The first is 1.5 times smaller
3. C.The first is 1.5 times larger
4. D.The first is 6 times larger
Show the answer
Correct answer: B. The first is 1.5 times smaller
$\eta = k/(\phi\mu c_t)$ : a factor of 2 in the numerator and 3 in the denominator give $2/3$ , so the hydraulic diffusivity of the first reservoir is 1.5 times smaller. Free gas slows down the transient processes more than the higher permeability speeds them up.
3
Estimate the relative error of the linear equation for oil with a compressibility of $2\cdot10^{-3}$ MPa⁻¹ at a pressure fall of 10 MPa, and for water with a compressibility of $0.45\cdot10^{-3}$ MPa⁻¹ at a fall of 20 MPa.
Show the answer
The error is about half of the product $c\,\Delta p$ . For oil $c\,\Delta p = 0.02$ and the error is about 1%. For water $c\,\Delta p = 0.009$ and the error is about 0.45%. In both cases the linear equation is more accurate than rate measurements.
4
A well is tested at a constant rate in a block that is bounded on all sides by sealing faults. Which conditions are set for the diffusivity equation?
1. A.Only a constant rate at the well
2. B.A constant pressure at the well and at the boundary
3. C.The initial pressure in the whole block, a constant rate at the well, zero flow across the boundary
4. D.A constant rate at the well and a constant pressure at the boundary
Show the answer
Correct answer: C. The initial pressure in the whole block, a constant rate at the well, zero flow across the boundary
Three conditions are needed: the initial one (the pressure equals $p_i$ everywhere), the one at the well (the pressure gradient that matches the rate), and the one at the outer boundary (through a sealing fault the flow is zero, $\partial p/\partial r = 0$ ). As long as the disturbance has not reached the faults, the solution coincides with the solution for an infinite reservoir. Then the pseudo-steady-state regime sets in.
5
Calculate the dimensionless time 10 h after the start of a well of radius 0.1 m, if $k$ = 20 mD, $\phi$ = 0.15, $\mu$ = 2 mPa·s and $c_t = 1.5\cdot10^{-3}$ MPa⁻¹.
Show the answer
$t_D = 3.553\cdot10^{-3}\,kt/(\phi\mu c_tr_w^2)$ . The numerator is $3.553\cdot10^{-3}\cdot20\cdot10 = 0.7106$ . The denominator is $0.15\cdot2\cdot1.5\cdot10^{-3}\cdot0.01 = 4.5\cdot10^{-6}$ . We get $t_D \approx 1.58\cdot10^{5}$ .
6
The permeability of a reservoir is 100 mD along the bedding and 10 mD across it. The axes of the reservoir simulation grid point along and across the bedding. Does such a model need a full permeability tensor with off-diagonal terms?
Show the answer
No. The axes coincide with the principal directions of permeability, so the tensor is diagonal: $k_x = k_y = 100$ mD, $k_z = 10$ mD, and the off-diagonal terms are zero. They would be needed if the grid were turned relative to the bedding, for example with cross-bedding.

Updated: 2 October 2026