Compressive VOF Module

Stream supports incompressible gas-liquid and multi-liquid interface capturing through the compressive_vof module. The module solves transport equations for liquid volume fractions (VOF), builds mixture density and viscosity from phase properties, and adds optional surface-tension and contact-angle forcing to the momentum equation.

Enable the module by adding loadModule: compressive_vof to the run control file:

loadModule: compressive_vof
{
  ... standard run control file content
}

Control File Setup

The module expects phase property definitions and VOF state definitions in the run control file.

Phase Properties

Two run-control variables define phase properties:

  • liquidPhases: one or more named incompressible liquid phases.

  • gasPhase: incompressible gas phase used for the non-liquid remainder.

Each phase must provide rho and mu.

Example phase setup
liquidPhases: <
  water=incompressible(rho=998.2 kg/m/m/m, mu=1.0e-3 kg/m/s),
  oil=incompressible(rho=850.0 kg/m/m/m, mu=4.0e-3 kg/m/s)
>
gasPhase: <
  air=incompressible(rho=1.225 kg/m/m/m, mu=1.8e-5 kg/m/s)
>

The liquid fraction for each phase is tracked separately. The gas fraction is computed as:

\[\alpha_{gas} = 1 - \sum_i \alpha_{liquid,i}\]

Initial Conditions

The module requires VOF initialization in initialCondition (uniform) or initialConditionRegions (regional states), using a vof=[phase=value,...] list whose phase names match those in liquidPhases.

Uniform VOF initial condition
initialCondition: <p=0 Pa, v=0.0 m/s, vof=[water=1.0, oil=0.0]>
Regional VOF initialization
initialConditionRegions: <
  default=state(vof=[water=0.0, oil=0.0]),
  slug=state(vof=[water=1.0, oil=0.0]),
  regions=[
    inBox(p1=[0.0, 0.0, 0.0], p2=[0.02, 0.002, 0.002], composition=slug)
  ]
>

Boundary Conditions

For incompressibleInlet, subsonicInlet, totalPressureInlet, inflow, and farfield boundaries, a vof list is required.

Inlet VOF specification
boundary_conditions: <
  Inlet=incompressibleInlet(v=[0.2 m/s, 0.0 m/s, 0.0 m/s], vof=[water=1.0, oil=0.0]),
  Outlet=fixedPressureOutlet(p=0.0 Pa),
  Wall=noslip
>

For no-slip walls, an optional per-boundary static contact angle override can be provided:

Wall=noslip(contactAngle=35 deg)

If contactAngle is provided, units must be degrees and the value must be in [0,180].

VOF Transport and Property Options

The following variables control VOF transport and property blending.

VOF and Mixture Property Controls

Variable

Allowed Values

Default

Description

vofInviscidFlux

FOU, IG, MIG, FBIC, TVD

FOU

Inviscid flux for liquid-volume-fraction transport.

vofViscosity

real

0.0

Extra numerical diffusion coefficient used in VOF assembly.

vofViscosityBlend

linear, harmonic

linear

Mixture laminar viscosity blending model.

overUnderFillingCorrection

on, off

off

Enables face-level over/under-fill correction in non-FOU VOF advection.

vofViscosityBlend=linear uses arithmetic blending by volume fraction. vofViscosityBlend=harmonic uses harmonic blending and is often preferable for high viscosity-ratio flows.

The solver also performs a final boundedness projection on each VOF solve collapse: individual liquid fractions are clipped to [0,1] and renormalized if their sum exceeds 1.

VOF Equation Solver Options

vofEquationOptions configures the VOF linear solve controls.

vofEquationOptions Entries

Entry

Allowed Values

Default

Description

linearSolver

SGS, FSGS, LSGS

SGS

Linear solver used for vofStar.

relaxationFactor

real in (0,1]

0.5

Relaxation factor for VOF equation update.

maxIterations

integer >= 0

5

Maximum VOF linear iterations per nonlinear iteration.

Example:

vofEquationOptions: <linearSolver=SGS, relaxationFactor=0.6, maxIterations=8>

Surface Tension, Curvature, and Contact Angle

Surface tension is activated by defining sigma.

Surface-Tension and Interface-Geometry Controls

Variable

Allowed Values

Default

Description

sigma

real (N/m)

not set

Surface tension coefficient.

vof_curvature_type

computed, none

computed

Curvature source: geometric computation or constant value.

kappaAnalytical

real

0.0

Constant curvature used when vof_curvature_type=none.

vof_curvature_regularization_type

none, clamp

none

Optional post-processing regularization of computed curvature.

vof_curvature_abs_max

real >0 (when clamp active)

1.0e6

Absolute curvature limit for clamp regularization.

vof_curvature_regularization_band

real in [0,0.5)

0.0

Apply clamp only for interface cells with band < capillaryColor < 1-band.

vof_contact_angle_model

none, static, dynamic

none

Wall contact-angle model for no-slip boundaries.

vof_contact_angle_deg

real in [0,180] deg

90.0

Baseline angle; for static model, this is the imposed value.

vof_dynamic_contact_angle_coeff

real >=0

9.0

Dynamic angle coefficient in Cox-Voinov-lite relation.

vof_dynamic_contact_angle_min_deg

real in [0,180] deg

5.0

Lower bound for dynamic contact angle.

vof_dynamic_contact_angle_max_deg

real in [0,180] deg

175.0

Upper bound for dynamic contact angle.

The dynamic contact-angle model uses:

\[\theta^3 = \theta_0^3 + C \cdot Ca\]

where \(\theta_0\) is from vof_contact_angle_deg, \(C\) is vof_dynamic_contact_angle_coeff, and \(Ca\) is a local wall-relative capillary number. The computed angle is then clamped by the min/max dynamic bounds.

If a no-slip boundary supplies contactAngle=..., that boundary value overrides the global vof_contact_angle_deg for that boundary patch. This override is active only when vof_contact_angle_model is static or dynamic.

Output and Diagnostics

When plot output is written, the module writes per-phase VOF scalar fields:

  • output/vof_<liquidName>_sca.<cycle>_<modelName>

  • output/vof_<gasName>_sca.<cycle>_<modelName>

Module scalar diagnostics:

  • output/liquidVolume.dat: total liquid volume for each liquid phase vs. cycle/time.

  • output/capillaryDiagnostics.dat with columns:

    • ncycle

    • stime

    • maxAbsKappa

    • maxCapillaryCFL

    • totalLiquidVolumeDriftAbs

    • totalLiquidVolumeDriftRel

The capillary CFL diagnostic is estimated as:

\[CFL_{\sigma} = \Delta t \sqrt{\frac{|\sigma|}{\rho V}}\]

evaluated near interface cells.

Complete Example

Minimal compressive_vof setup for an air-water case
loadModule: compressive_vof
{
  liquidPhases: <water=incompressible(rho=998.2 kg/m/m/m, mu=1.0e-3 kg/m/s)>
  gasPhase: <air=incompressible(rho=1.225 kg/m/m/m, mu=1.8e-5 kg/m/s)>

  initialCondition: <p=0 Pa, v=0.0 m/s, vof=[water=0.0]>

  boundary_conditions: <
    Inlet=incompressibleInlet(v=[0.1 m/s, 0.0 m/s, 0.0 m/s], vof=[water=1.0]),
    TubeWall=noslip(contactAngle=35 deg),
    Outlet=fixedPressureOutlet(p=0.0 Pa)
  >

  vofEquationOptions: <linearSolver=SGS, relaxationFactor=0.6, maxIterations=8>
  vofInviscidFlux: MIG
  vofViscosityBlend: harmonic
  overUnderFillingCorrection: on

  sigma: 0.072 N/m
  vof_curvature_type: computed
  vof_curvature_regularization_type: clamp
  vof_curvature_abs_max: 2.0e4
  vof_contact_angle_model: dynamic
  vof_contact_angle_deg: 90.0
  vof_dynamic_contact_angle_coeff: 9.0
  vof_dynamic_contact_angle_min_deg: 20.0
  vof_dynamic_contact_angle_max_deg: 160.0
}

Notes

  • The current surface-tension model uses a single sigma value for all interfaces.

  • Multiple liquid phases are supported in transport and property blending.

  • For gasPhase, if multiple gas entries are provided, the last one read is used for rhoGas and muGas.