Multiphase, Interfacial, and Free-Surface Flows

Introduction and Context

Multiphase flow occurs when two or more physically distinct phases—such as liquid–gas, liquid–solid, or liquid–liquid—move simultaneously through a domain. Examples include boiling in steam generators, cavitation in hydraulic machinery, sediment transport in rivers, oil–water flow in pipelines, pneumatic conveying, spray combustion, and blood flow containing suspended cells. Free-surface flow is a special class in which the interface between a liquid and a gas is deformable and exposed to a known pressure, commonly atmospheric pressure.

The defining difficulty is that the interface is not merely a boundary: it moves, deforms, exchanges mass and momentum, and may undergo topological changes such as breakup, coalescence, entrainment, and phase change. Consequently, multiphase-flow analysis combines conservation laws with interfacial thermodynamics, constitutive modelling, dimensional analysis, and numerical interface tracking.

At postgraduate level, it is essential to distinguish between a mathematically exact description and an engineering closure model. The Navier–Stokes equations can describe each phase, but practical prediction also requires models for turbulence, interfacial momentum transfer, wall forces, phase change, droplet or bubble size distributions, and unresolved interface geometry.

Learning Objectives

By the end of this chapter, the reader should be able to:

1. Classify multiphase and free-surface flows using physically meaningful criteria.
2. Formulate conservation equations for phases and interfaces.
3. Explain surface tension, capillarity, wetting, and contact-line motion.
4. Use dimensionless groups to identify dominant mechanisms.
5. Interpret flow regimes in channels, pipes, and particle-laden systems.
6. Select and critically evaluate interface-capturing, interface-tracking, and averaged models.
7. Analyse a practical free-surface or multiphase-flow problem and identify modelling limitations.

1. Terminology and Classification

A phase is a region of matter with relatively uniform composition and thermodynamic state. A component is a chemical constituent. Thus, water and oil are two phases and usually two components, whereas liquid water and water vapour are two phases of one component.

A dispersed phase exists as bubbles, droplets, or particles within a continuous phase. In a bubbly liquid, gas is dispersed and liquid is continuous. In a spray, liquid droplets are dispersed in gas. In stratified pipe flow, both phases may be continuous over substantial regions.

The volume fraction of phase k is

where V_k is the volume occupied by phase k inside a representative volume V. For phases that completely fill the control volume,

The void fraction is the gas volume fraction in a gas–liquid system. The superficial velocity of phase k is

where Q_k is the volumetric flow rate that phase would have if it alone occupied the full cross-sectional area A. The phase-average velocity differs from superficial velocity:

provided α_k is nonzero.

Multiphase flows may be classified by:

• phase combination: gas–liquid, liquid–solid, gas–solid, or liquid–liquid;
• spatial arrangement: dispersed, stratified, annular, slug, or interpenetrating;
• flow regime: laminar, transitional, or turbulent;
• interface behaviour: smooth, wavy, breaking, or topologically changing;
• thermodynamic behaviour: non-evaporating, evaporating, condensing, reactive, or cavitating.

The classification is not merely descriptive. It determines which averaged variables, interfacial closures, and numerical methods are credible.

2. Governing Conservation Laws

2.1 Single-phase local equations

For an incompressible Newtonian fluid of density ρ and dynamic viscosity μ, the continuity and momentum equations are

Here u is velocity, p is pressure, g is gravitational acceleration, and f_s represents body-force density or a singularly distributed interfacial force. The equations assume a continuum description and constant density; the momentum equation can be generalized to compressible or non-Newtonian phases.

2.2 Two-fluid formulation

In an Eulerian two-fluid model, each phase is treated as an interpenetrating continuum. For phase k,

ρ_k, u_k, p_k, and τ_k are phase density, velocity, pressure, and viscous stress tensor. Γ_k is the mass-transfer rate into phase k per unit volume, M_k is interfacial momentum transfer, and u_Γ is the velocity associated with transferred mass. Conservation requires

and, in the absence of external interphase forces,

The two-fluid model is powerful for engineering systems but requires closure laws for drag, lift, added mass, turbulent dispersion, wall lubrication, and interfacial area density.

2.3 Interface jump conditions

Let n point from phase 1 to phase 2 and let [q] = q_2 − q_1 denote the jump in a quantity q. For an interface moving with normal velocity V_n, mass conservation gives

For a non-evaporating interface, the normal velocity is continuous in the appropriate material sense. The stress balance is

where T is the Cauchy stress tensor, σ is surface tension, κ is the signed curvature, and ∇_s is the surface gradient. For constant σ, the tangential Marangoni term ∇_sσ vanishes. The normal component yields the Young–Laplace pressure jump:

For a spherical interface of radius R, κ = 2/R, so

This result assumes a static spherical interface, negligible viscous stress, and uniform surface tension. For a soap film with two interfaces, the pressure difference is approximately 4σ/R.

3. Surface Tension, Wetting, and Capillarity

Surface tension is the excess free energy per unit interfacial area and has units N/m, equivalent to J/m². It arises because molecules at an interface experience an asymmetric molecular environment. Surface tension tends to reduce area, which explains why small droplets and bubbles become nearly spherical.

At a solid–liquid–gas contact line, the equilibrium contact angle θ_e is related to interfacial energies by Young’s equation:

The symbols σ_SG, σ_SL, and σ_LG denote solid–gas, solid–liquid, and liquid–gas interfacial tensions. This equation assumes an ideal smooth, chemically homogeneous, and stationary surface. Real surfaces exhibit roughness, chemical heterogeneity, contact-angle hysteresis, pinning, and dynamic contact angles.

A useful capillary length is

where ∆ρ is the relevant density difference. Disturbances much smaller than ℓ_c are strongly influenced by surface tension; larger disturbances are more strongly influenced by gravity.

For a capillary tube of radius r, the equilibrium height change for a liquid of contact angle θ is

This expression assumes a circular tube, hydrostatic equilibrium, negligible evaporation, and a static contact angle. It shows why capillary effects become important in porous media, microfluidics, heat pipes, and biological flows.

4. Dimensionless Groups and Physical Regimes

Dimensionless groups compare competing mechanisms.

The Reynolds number compares inertia and viscosity:

where U and L are characteristic velocity and length, and ν = μ/ρ is kinematic viscosity.

The Weber number compares inertia and surface tension:

The Froude number compares inertia and gravity:

The Bond number compares gravity and surface tension:

For a density difference, ∆ρ replaces ρ when buoyancy is the relevant force. The capillary number compares viscous stress and surface tension:

The Ohnesorge number combines viscosity, inertia, and surface tension:

A low Weber number indicates that surface tension resists deformation; a high Weber number favours stretching and breakup. A low Froude number indicates gravity-dominated free-surface motion, while a high Froude number indicates inertia-dominated motion. These interpretations require a consistent choice of U and L; different local scales can produce different regime classifications.

5. Free-Surface Flows

A free surface is an interface whose pressure is prescribed, usually p_atm. For a clean, non-evaporating liquid exposed to gas, the normal-stress condition is approximately

If gas density and viscosity are negligible and curvature is small, p_liquid is approximately atmospheric at the surface.

The kinematic boundary condition states that fluid particles on the interface remain on it. If the interface is represented by F(x,y,z,t) = 0, then

for a material interface. If F is a level-set function, its zero contour is the interface. The dynamic condition specifies stress balance, while the kinematic condition specifies interface motion.

Shallow-water approximation

For a long, shallow channel with depth h much smaller than wavelength L, vertical acceleration is small and pressure is approximately hydrostatic. The one-dimensional Saint-Venant equations for width b constant are

Here u is depth-averaged velocity, S_0 is bed slope, and S_f is friction slope. The equations neglect vertical structure, dispersive effects, and strong curvature. They are widely used for floods, open-channel transients, hydraulic jumps, and dam-break approximation.

The Froude number for open-channel flow is

Subcritical flow has Fr < 1, critical flow has Fr = 1, and supercritical flow has Fr > 1. Information propagates differently in these regimes, which affects boundary-condition specification.

Hydraulic jump

A hydraulic jump converts supercritical flow into subcritical flow with substantial energy dissipation. For a rectangular channel, the sequent-depth relation is

where y₁ is the upstream depth, y₂ is the downstream depth, and Fr₁ is the upstream Froude number. The relation assumes steady, one-dimensional flow, hydrostatic pressure distribution away from the jump, negligible bed slope over the jump length, and conservation of momentum rather than mechanical energy.

6. Multiphase Flow Regimes

Gas–liquid systems

Common pipe-flow patterns include bubbly flow, slug flow, churn flow, annular flow, and stratified flow. Regime transitions depend on superficial velocities, density ratio, viscosity ratio, pipe diameter, inclination, and surface tension. Flow maps are therefore empirical or semi-empirical and should not be transferred uncritically between fluids or geometries.

In dispersed bubbly flow, interfacial drag is often represented through a drag coefficient C_D. For a spherical particle or bubble of diameter d moving with slip speed U_r, the drag magnitude is commonly written

where ρ_c is continuous-phase density and A_p = πd²/4 is projected area. The value of C_D depends on Reynolds number, shape, contamination, deformation, and proximity to walls or other bubbles. Using an isolated rigid-sphere correlation for deformable bubbles can produce serious errors.

Liquid–solid systems

Particle settling results from the balance of gravity, buoyancy, and drag. For a small spherical particle in creeping flow, the Stokes settling velocity is

where ρ_p and ρ_f are particle and fluid densities, d is particle diameter, and μ_f is fluid viscosity. The assumptions are spherical particles, very small particle Reynolds number, dilute suspension, and an unbounded Newtonian fluid. The associated particle Reynolds number is

If Re_p is not much less than unity, inertial drag corrections are required.

Liquid–liquid systems

Liquid–liquid flow may be stratified, dispersed, or emulsified. Interfacial tension, viscosity ratio, density ratio, surfactants, and residence time control droplet deformation and coalescence. Surfactants may lower σ but can also generate Marangoni stresses that resist interface drainage and alter coalescence kinetics.

7. Modelling and Numerical Methods

Interface tracking

In front-tracking methods, the interface is represented explicitly by connected marker points or surface elements. These methods provide accurate geometric information and curvature but require special treatment when interfaces merge or break.

Interface capturing

The volume-of-fluid (VOF) method transports a phase volume fraction:

for incompressible advection without phase change. The interface is reconstructed from cell fractions. VOF conserves volume well, which is valuable for long simulations, but curvature and surface-tension calculation can suffer from parasitic currents if the pressure and curvature discretizations are inconsistent.

The level-set method transports a signed-distance-like function φ:

The interface is φ = 0. Level sets provide smooth normals and curvature but may lose or gain mass unless a correction or conservative formulation is used.

The phase-field method introduces a diffuse order parameter and solves an additional free-energy-based evolution equation. It naturally handles topological change and wetting, but the interface thickness must be resolved and its parameters require physical interpretation.

Averaged and Eulerian–Lagrangian approaches

Eulerian–Eulerian models solve averaged equations for all phases and are suitable for large industrial domains. Eulerian–Lagrangian models solve the carrier phase on a grid and track particles or droplets individually. The latter is effective for dilute dispersed flows but becomes expensive or statistically unreliable at high volume fraction.

A model should be selected according to the required output. If the objective is pressure drop in a large pipe, a validated averaged model may be appropriate. If the objective is ligament breakup or contact-line dynamics, an interface-resolving method is usually necessary.

8. Turbulence, Closure, and Uncertainty

Averaging nonlinear equations generates correlations such as Reynolds stresses, turbulent scalar fluxes, and fluctuating interfacial forces. Multiphase turbulence is not obtained simply by applying a single-phase turbulence model independently to each phase. Interphase slip, wakes, bubble-induced turbulence, preferential concentration, and modulation of turbulence must be considered.

Important closure quantities include:

• interfacial area concentration;
• drag, lift, and virtual-mass coefficients;
• wall-film thickness and entrainment rate;
• collision and coalescence kernels;
• breakup frequency and daughter-size distribution;
• phase-change mass flux;
• dynamic contact-angle laws.

Validation should compare quantities at the scale of intended use. Agreement in mean pressure drop does not establish accuracy in instantaneous interface topology. Grid convergence, sensitivity to closure parameters, conservation checks, and uncertainty quantification are essential parts of a credible simulation.

9. Worked Examples

Example 1: Capillary rise

Water at room temperature has σ = 0.072 N/m, ρ = 998 kg/m³, and contact angle θ = 20° in a clean capillary tube of radius r = 0.50 mm. Estimate the equilibrium rise.

Using

and cos 20° ≈ 0.940,

Thus, h ≈ 27.7 mm. The result neglects evaporation, contamination, meniscus hysteresis, and dynamic effects. The units are consistent because N/m divided by N/m² gives metres.

Example 2: Hydraulic-jump depth

Water flows in a rectangular channel with upstream depth y₁ = 0.20 m and upstream Froude number Fr₁ = 3.0. Determine the sequent depth.

Therefore,

Hence,

The large increase in depth is consistent with a strong hydraulic jump. The equation predicts momentum-equivalent depth, not a loss-free energy transformation.

10. Applied Case: Slugging in a Gas–Liquid Pipeline

Background

Consider an illustrative offshore pipeline carrying natural gas and condensate from a subsea manifold to a separator. The line is inclined, contains low points, and experiences fluctuating production rates. Operators observe periodic surges in separator liquid level and pressure.

Problem

The likely mechanism is severe slugging: liquid accumulates at a low point, blocks gas passage, and is periodically expelled as a large liquid slug. The event creates high instantaneous liquid flow, separator overload, vibration, and possible control-valve instability.

Analysis

A simple homogeneous model using average mixture density may predict mean pressure drop but miss the intermittent accumulation. A regime-oriented analysis considers superficial gas and liquid velocities, inclination, hydrostatic head, compressibility of the gas pocket, and outlet backpressure. Low gas velocity and strong inclination favour liquid holdup and slug formation. A transient multiphase model is required because the system stores mass and compressible gas energy.

Approach

A practical engineering workflow is:

1. Establish fluid properties over the operating pressure and temperature range.
2. Construct a transient model with segmented pipe geometry.
3. Compare predicted slug frequency, liquid volume, pressure amplitude, and arrival time with field data.
4. Test operating changes such as choke adjustment, gas-lift injection, inlet-flow smoothing, or a slug catcher.
5. Verify that the proposed control strategy remains stable under start-up and shut-down conditions.

A slug catcher or adequately sized separator reduces the consequence but does not necessarily eliminate the hydrodynamic cause. Increasing gas velocity may suppress severe slugging in some systems, but it can increase erosion, compressor demand, and pressure loss.

Outcome

For this illustrative case, a combined strategy—controlled inlet choking, improved level-control tuning, and installation of a slug catcher—reduces separator excursions while preserving acceptable production. The model is accepted only after reproducing measured transients within an agreed engineering tolerance.

Learning points

The case demonstrates that a steady pressure-drop calculation is insufficient for an unsteady multiphase system. Geometry, storage, compressibility, closure models, and control dynamics must be considered together. It also illustrates why validation data are essential: different models may predict similar mean conditions but very different transient risks.

11. Practical Design and Best Practices

A robust multiphase-flow investigation should:

1. Define the required prediction: mean pressure drop, interface location, transient load, heat-transfer coefficient, or failure probability.
2. Select characteristic scales and calculate relevant dimensionless groups.
3. Determine whether phases are continuous or dispersed and whether topology changes matter.
4. Use conservation-preserving numerical schemes.
5. Resolve the interface and smallest dynamically important scales where feasible.
6. Use closure correlations within their validated ranges.
7. Perform mesh, time-step, and parameter sensitivity studies.
8. Validate against experiments or plant data at operating conditions resembling the application.
9. Report uncertainty, not only a single deterministic prediction.
10. Check erosion, corrosion, cavitation, pressure relief, structural vibration, and environmental consequences.

12. Limitations, Safety, and Contemporary Developments

Multiphase predictions are often less certain than single-phase predictions because interfaces are unresolved, measurements are difficult, and closure correlations are system-dependent. Contact-line physics, turbulent breakup, surfactant transport, and phase-change kinetics remain active research areas.

High-fidelity simulation is increasingly combined with adaptive meshing, reduced-order models, experimental data assimilation, and machine-learning surrogates. Artificial intelligence can accelerate parameter estimation or regime classification, but it does not replace conservation laws, dimensional reasoning, or validation. A data-driven model trained on one fluid, geometry, or scale may fail under extrapolation.

Safety considerations include overpressure from blocked flow, water hammer, cavitation damage, vessel rollover, toxic releases, erosion by particles, and loss of containment. Designs must comply with applicable pressure-vessel, process-safety, environmental, and occupational regulations. Ethical engineering practice requires transparent reporting of assumptions and uncertainty, particularly where model outputs influence safety margins.

Chapter Summary

Multiphase and free-surface flows involve interacting phases, moving interfaces, and often unresolved geometric and thermodynamic processes. The volume fraction and superficial velocity provide useful phase descriptions, while interface jump conditions express mass and momentum conservation at a deformable boundary. Surface tension produces the Young–Laplace pressure jump, controls capillary length, and interacts with wetting and contact-angle physics.

Reynolds, Weber, Froude, Bond, capillary, and Ohnesorge numbers identify the relative importance of inertia, viscosity, gravity, and surface tension. Free-surface flows can often be simplified through hydrostatic or shallow-water assumptions, but those assumptions fail near breaking waves, hydraulic jumps, strong curvature, and rapid transients. Multiphase pipe and particle-laden flows require regime classification and closure laws for interphase forces and population dynamics.

VOF, level-set, front-tracking, phase-field, Eulerian–Eulerian, and Eulerian–Lagrangian methods each involve trade-offs in conservation, geometry, computational cost, and treatment of topology. Reliable engineering analysis requires appropriate model selection, validation, sensitivity studies, and explicit uncertainty assessment. The pipeline slugging case illustrates that transient storage, compressibility, geometry, and control must be analysed together.

Key Terms

Phase: A physically distinct region with relatively uniform composition and state.

Multiphase flow: Flow involving two or more phases.

Free surface: A deformable liquid–gas interface exposed to a prescribed gas pressure.

Volume fraction: Fraction of a representative volume occupied by a specified phase.

Void fraction: Gas volume fraction in a gas–liquid system.

Superficial velocity: Phase flow rate divided by total cross-sectional area.

Slip velocity: Relative velocity between two phases.

Surface tension: Interfacial free energy per unit area, measured in N/m.

Contact angle: Angle characterizing liquid wetting at a solid surface.

Young–Laplace equation: Relation between pressure jump and interface curvature.

Capillary length: Length scale at which gravity and surface tension are comparable.

Reynolds number: Ratio of inertial to viscous effects.

Weber number: Ratio of inertial to surface-tension effects.

Froude number: Ratio of inertial to gravitational effects.

Bond number: Ratio of gravitational to surface-tension effects.

VOF method: Conservative interface-capturing method based on volume fraction.

Level-set method: Interface method using a transported scalar whose zero contour defines the interface.

Interfacial area concentration: Interface area per unit mixture volume.

Hydraulic jump: Abrupt transition from supercritical to subcritical open-channel flow.

Severe slugging: Large-amplitude intermittent gas–liquid flow caused by accumulation and release of liquid.

Review Questions

Short Answer Questions

1. Distinguish between volume fraction, phase-average velocity, and superficial velocity.
2. State the assumptions behind the Young–Laplace pressure relation for a spherical interface.
3. What physical mechanisms are compared by the Weber and Froude numbers?
4. Why can a two-fluid model not be used without interfacial closure relations?
5. State one principal advantage and one limitation of the VOF method.

Long Answer Questions

1. Derive and explain the kinematic and dynamic boundary conditions for a clean, non-evaporating free surface.
2. Compare front-tracking, VOF, level-set, and phase-field methods in terms of conservation, curvature calculation, and topological change.
3. Explain the principal gas–liquid flow regimes in a pipe and discuss the variables controlling regime transition.
4. Derive the shallow-water continuity and momentum equations and state their assumptions.
5. Discuss how turbulence and unresolved interfacial geometry create closure problems in averaged multiphase models.

Analytical / Application Questions

1. A capillary tube has radius 0.25 mm, water surface tension 0.072 N/m, density 998 kg/m³, and contact angle 30°. Calculate the equilibrium capillary rise and state the main assumptions.
2. A rectangular channel has upstream depth 0.15 m and upstream Froude number 2.5. Calculate the sequent depth after a hydraulic jump and interpret the result.
3. An oil–water pipeline experiences periodic pressure surges and separator-level oscillations. Propose a diagnostic and modelling workflow, identifying measurements, relevant dimensionless groups, possible flow regimes, and validation criteria.

Multiple Choice Questions with Answers

Question 1

Which quantity is defined as the phase flow rate divided by the total pipe cross-sectional area?

A. Phase-average velocity

B. Superficial velocity

C. Slip velocity

D. Terminal velocity

Correct Answer: B. Superficial velocity

Question 2

For a static spherical interface with uniform surface tension, the pressure jump is given by which relation?

A. ∆p = σR/2

B. ∆p = 2σ/R

C. ∆p = ρgR

D. ∆p = μU/R

Correct Answer: B. ∆p = 2σ/R

Question 3

The Weber number primarily compares inertia with:

A. Viscosity

B. Gravity

C. Surface tension

D. Compressibility

Correct Answer: C. Surface tension

Question 4

Which condition represents supercritical open-channel flow?

A. Fr < 1

B. Fr = 0

C. Fr = 1

D. Fr > 1

Correct Answer: D. Fr > 1

Question 5

For a non-evaporating material interface represented by F = 0, which equation expresses the kinematic condition?

A. ∂F/∂t + u · ∇F = 0

B. ∇ · u = 1

C. p = ρg

D. ∂p/∂t = 0

Correct Answer: A. ∂F/∂t + u · ∇F = 0

Question 6

Which interface method is directly based on transporting cell phase fractions?

A. Front tracking

B. Volume of fluid

C. Boundary element method

D. Particle image velocimetry

Correct Answer: B. Volume of fluid

Question 7

The Stokes settling-velocity formula is valid primarily when the particle Reynolds number is:

A. Much less than one

B. Approximately one hundred

C. Equal to the Weber number

D. Much greater than one thousand

Correct Answer: A. Much less than one

Question 8

A hydraulic jump changes the flow from:

A. Subcritical to supercritical

B. Laminar to turbulent only

C. Supercritical to subcritical

D. Compressible to incompressible

Correct Answer: C. Supercritical to subcritical

Question 9

Which dimensionless number compares viscous stress with surface tension?

A. Capillary number

B. Froude number

C. Bond number

D. Mach number

Correct Answer: A. Capillary number

Question 10

A major limitation of a steady homogeneous pressure-drop model for severe slugging is that it neglects:

A. Fluid density

B. Mass storage and transient compressibility effects

C. Pipe diameter

D. Gravity entirely in every formulation

Correct Answer: B. Mass storage and transient compressibility effects