Finite Volume Method — HEM

The wellbore simulator solves the one-dimensional area-weighted conservation laws for transient CO₂ flow using a cell-centred finite volume method. This page describes the Homogeneous Equilibrium Model (HEM) discretisation, following Sacconi & Mahgerefteh (2020). Two-phase regions are treated as a single pseudo-fluid in mechanical and thermodynamic equilibrium (one velocity, one pressure, one temperature).

Governing Equations

For variable cross-section $A(z)$ and inclination $\theta$ (angle from horizontal, so a vertical downward well has $\theta = -90°$), the area-weighted conservation laws read:

$$\frac{\partial (\rho A)}{\partial t} + \frac{\partial (\rho u A)}{\partial z} = 0 \tag{1}$$ $$\frac{\partial (\rho u A)}{\partial t} + \frac{\partial (\rho u^2 A + p A)}{\partial z} = p \, \frac{\partial A}{\partial z} + A\bigl(F - \rho g \sin\theta\bigr) \tag{2}$$ $$\frac{\partial (\rho E A)}{\partial t} + \frac{\partial \bigl[(\rho E + p) u A\bigr]}{\partial z} = A\bigl(F u - \rho u g \sin\theta + Q\bigr) \tag{3}$$

with total energy $E = e + \tfrac{1}{2}u^2$, friction force per unit volume $F$, and wall heat source per unit volume $Q$ (positive into fluid).

State Vector

The conserved state stored per cell is area-weighted:

$$\mathbf{U} = \begin{bmatrix} \rho A \\ \rho u A \\ \rho E A \\ A \end{bmatrix} \tag{4}$$

Cell-Centred Finite Volume Update

Discretising the wellbore into cells of size $\Delta z_i$, integrating (1)–(3) over cell $i$, and applying the divergence theorem yields the semi-discrete update:

$$\frac{d\mathbf{U}_i}{dt} = -\frac{\mathbf{F}_{i+\frac{1}{2}} - \mathbf{F}_{i-\frac{1}{2}}}{\Delta z_i} + \mathbf{S}_i^{\text{geom}} + \mathbf{S}_i \tag{5}$$

with three contributions on the right-hand side:

• $\mathbf{F}_{i\pm\frac{1}{2}}$ — numerical fluxes at cell faces from the area-weighted AUSM⁺-up scheme. See the AUSM⁺-up Flux Scheme page for splitting polynomials, well-balanced face pressures, and dissipation terms. Higher-order accuracy comes from MUSCL or WENO5 reconstruction of left/right interface states.
• $\mathbf{S}_i^{\text{geom}}$ — non-conservative geometric pressure work, contributing $p \, \partial A/\partial z$ to the momentum equation. Evaluated from face pressures consistent with the AUSM⁺-up reconstruction.
• $\mathbf{S}_i$ — physical source terms (gravity, friction, heat transfer), assembled below.

Source Terms

In HEM the source vector is

$$\mathbf{S}_i = \begin{bmatrix} 0 \\ A_i\,(F_i - \rho_i g \sin\theta_i) \\ A_i\,(F_i u_i - \rho_i u_i g \sin\theta_i + Q_i) \\ 0 \end{bmatrix} \tag{6}$$

with the closures:

Friction (Darcy–Weisbach with Chen's correlation):

$$F = -\frac{f_D\, \rho\, u\,|u|}{2 D}, \qquad f_D = f_D^{\text{Chen}}\bigl(\mathrm{Re}, \varepsilon/D\bigr) \tag{7}$$

where $f_D$ is the Darcy friction factor from Chen (1979), evaluated from the Reynolds number $\mathrm{Re} = \rho |u| D / \mu$ and relative roughness $\varepsilon/D$. Viscosity $\mu$ comes from the equation of state.

Wall heat source: $Q$ is computed by a transient radial conduction model coupling fluid → tubing → annulus → casing → cement → formation, evaluated per cell from the local wall and formation temperatures.

Both gravity and friction contributions are individually toggleable via input flags.

Closure

The HEM model requires a thermodynamic closure $p = p(\rho, e)$ valid across single-phase liquid, single-phase vapour, and the two-phase saturation dome. The simulator obtains this from CoolProp, evaluated either directly per call or via pre-tabulated lookup for performance. Inside the saturation envelope, pressure is fixed at $p_{\mathrm{sat}}(T)$ and the mixture density follows from the lever rule on the conserved internal energy.

Time Integration

Time integration is delegated to PETSc TS, with a Newton solver (SNES) for implicit schemes. Step size is adapted from a CFL bound using the maximum interface wave speed $|u| + a$ returned by the AUSM⁺-up flux routine.

References

• Sacconi, A. & Mahgerefteh, H. (2020). Modelling start-up injection of CO₂ into highly-depleted gas fields. Energy, 191, 116530. doi:10.1016/j.energy.2019.116530
• Chen, N. H. (1979). An explicit equation for friction factor in pipe. Industrial & Engineering Chemistry Fundamentals, 18(3):296–297.
• Liou, M.-S. (2006). A sequel to AUSM, Part II: AUSM⁺-up for all speeds. Journal of Computational Physics, 214(1):137–170. doi:10.1016/j.jcp.2005.09.020