Conjugate Scalar Transport

• Conjugate scalar transport is the modeling of scalar fields across fluid and solid phases with prescribed discontinuities at sharp interfaces.
• It employs a volume penalization method with diffuse-interface representations to seamlessly integrate advection–diffusion in fluids and pure diffusion in solids.
• Validation studies show that the method captures both scalar and flux jumps with relative errors below 3%, ensuring reliable simulations in complex geometries.

Conjugate scalar transport governs the coupled evolution of a scalar field—such as temperature or concentration—in multiphase domains, typically involving both fluid and solid regions separated by a sharp interface. This phenomenon is central to a wide range of thermal and chemical processes wherein the scalar and its normal flux may exhibit prescribed discontinuities (“jump conditions”) across the interface. Mathematical and computational treatments of conjugate scalar transport must systematically unify advection–diffusion equations with diverse interfacial boundary conditions over complex geometries (Liu et al., 15 Jan 2026).

1. Governing Equations for Multiphase Domains

Given a domain $\Omega_e = \Omega_f \cup \Omega_s$ partitioned into fluid ($\Omega_f$) and solid ($\Omega_s$), distinct scalar fields $c_f$ and $c_s$ are defined for each phase. The canonical governing equations are:

• In the fluid (advection–diffusion):

$$\partial_t c_f + \mathbf{u} \cdot \nabla c_f = \nabla \cdot (D_f \nabla c_f)$$

• In the solid (pure diffusion):

$\partial_t c_s = \nabla \cdot (D_s \nabla c_s)$

where $\mathbf{u}$ is the velocity, and $D_f$, $D_s$ are the molecular diffusivities in each region.

At the sharp interface $\Gamma = \partial\Omega_f \cap \partial\Omega_s$, the following generalized interfacial boundary conditions are imposed:

$[c] \equiv c_s - c_f = g_1(x, t)$

$$[D\nabla c \cdot \mathbf{n}] \equiv D_s \nabla c_s \cdot \mathbf{n} - D_f \nabla c_f \cdot \mathbf{n} = g_2(x, t)$$

where $\mathbf{n}$ is the outward unit normal from solid to fluid and $g_1$, $g_2$ are prescribed jumps in scalar and normal flux, respectively.

2. Diffuse-Interface Representation and Penalization Methodology

Simulations of conjugate scalar transport involving arbitrary interface geometry are computationally challenging, particularly when accommodating jump conditions. The volume penalization method (VPM), an immersed boundary approach, provides a framework for unifying the treatment of disparate subdomains and enforcing interfacial jump conditions (Liu et al., 15 Jan 2026).

A signed-distance level-set function $\phi_0(x)$ encodes the interface: $\phi_0 = 0$ on $\Gamma$, $\phi_0 > 0$ in $\Omega_s$, $\phi_0 < 0$ in $\Omega_f$. A smoothed phase indicator $\phi(x)$ transitions from $0$ (fluid) to $1$ (solid) over a thickness $2\delta$, defined as

$$\phi(x) = \begin{cases} 0, & \phi_0 < -\delta \ \frac{1}{2}\left[1 + \sin\left(\frac{\pi \phi_0}{2\delta}\right)\right], & |\phi_0| \leq \delta \ 1, & \phi_0 > +\delta \end{cases}$$

The unit normal is given by $\mathbf{n} = -\nabla \phi_0 / |\nabla \phi_0|$.

For problems with a flux jump $g_2 = q_w$ (pure Neumann), the penalization introduces a source confined to the interfacial region:

$$S_{\text{jump}}(x) = \nabla \cdot [q_w\, \phi(x)\, \mathbf{n}(x)]$$

This source asymptotically approaches a Dirac delta at $\Gamma$ as $\delta \to 0$, enforcing the required flux jump in the integral sense.

3. Unified Penalized Equation and Extension to Scalar Jumps

The volume penalization approach enables the solution of a unified equation across $\Omega_f \cup \Omega_s$, with space-dependent coefficients:

$$\partial_t c + \mathbf{u} \cdot \nabla c = \nabla \cdot [D(x) \nabla c] + S_{\text{jump}}(x)$$

where

$$D(x) = \begin{cases} D_f, & \phi_0 < 0 \ D_s, & \phi_0 > 0 \end{cases}$$

If both scalar and flux jumps exist ($g_1 \neq 0$, $g_2 \neq 0$), an equivalent scalar field $h(x, t) = \alpha(x) c(x, t)$ is defined, with

$$\alpha(x) = \begin{cases} 1, & \phi_0 < 0 \ \alpha_s, & \phi_0 > 0 \end{cases}$$

such that $[h] = 0$ on $\Gamma$. The effective interfacial flux jump is

$$q_w^* = D(x)\left(1 - \frac{1}{\alpha(x)}\right) \partial_n h + g_2(x, t)$$

and the final penalized system becomes

$$\partial_t h + \mathbf{u} \cdot \nabla h = \nabla \cdot [D(x)\nabla h] + \nabla \cdot [q_w^*(x)\, \phi(x)\, \mathbf{n}(x)]$$

Recovery of $c$ is achieved via $c = h/\alpha$ (Liu et al., 15 Jan 2026).

4. Discretization and Computational Implementation

The VPM for conjugate scalar transport utilizes a finite-volume discretization on a Cartesian mesh. Second-order central differences are employed for both advection and diffusion operators. Time integration is handled by the fully implicit Euler scheme (first order in time). The interface mask, including the phase indicator and unit normal, is computed from the analytic or reinitialized level-set $\phi_0$; typically, $\delta = K_\delta \Delta$ with $K_\delta \approx 1.5$ (where $\Delta$ is the grid spacing). Adaptive mesh refinement (AMR) can be activated in regions where $|\phi_0| < 2\delta$ to ensure an adequate number of cells across the interface.

Outline of the VPM Algorithm

Given signed-distance φ₀(x) on the Cartesian mesh: 1. Compute phase-indicator φ(x) and unit normal n(x) = −∇φ₀/|∇φ₀|. 2. Set D(x), α(x) from φ₀. 3. Initialize h(x,0) = α(x)·c₀(x). Loop over time steps: a. Recompute ∂ₙh ≈ n · ∇h. b. Compute interfacial flux jump source: S_jump = ∇·[q_w*(x) · φ(x) · n(x)], where q_w*(x) = D(x)(1−1/α(x))·∂ₙh + prescribed g₂(x,t). c. Assemble and solve (implicit Euler): (h^{n+1} − h^n)/Δt + u·∇h^{n+1} = ∇·[D ∇h^{n+1}] + S_jump^{n+1}. d. Recover c^{n+1} = h^{n+1}/α. end loop

5. Validation, Error Analysis, and Comparison

The VPM has been verified in several canonical test problems:

• 1D diffusion with pure Neumann boundary: For $y' \in [0, d]$, analytical solution $c = 2 D_f t + y'^2$ is recovered with first-order grid convergence. Average relative error is approximately $1.2\%$ for $\Delta x = 0.02$, improving over the previous method (Brown-Dymkoski et al. 2014) which incurred up to $3\%$ error and exhibited nonphysical $c > 0$ in the solid.
• 2D fluid–solid coupled diffusion: For a cylindrical geometry ($r=0.1$ in $[-0.5,0.5]^2$), four cases with varying $(D_f, D_s, \alpha_s, q_w)$ closely match reference body-fitted mesh (BFM) simulations; all relative deviations along $y=0$ are below $3\%$ (mean $1.46\%$). Mesh convergence on Cartesian grids with AMR achieves less than $0.5\%$ change for $\geq 200 \times 200$ cells.
• 2D advection–diffusion: At $Re=40$, $Sc=0.5$, $\kappa = D_s/D_f = 5$, $\alpha_s=1.25$, $q_w=100$ on a cylinder of $r=0.5$ in $[-5,10]\times[-5,5]$, the VPM achieves a relative deviation of approximately $1.59\%$ versus BFM (chtMultiRegionFoam). Both scalar and flux jumps are captured.

In all test cases, the VPM enforces the interfacial jump in an integral sense while confining the penalization source to the narrow interfacial band, thereby avoiding pollution of the solid region (Liu et al., 15 Jan 2026).

6. Advantages, Limitations, and Interpretations

The VPM presented provides a unified diffuse-interface formulation for fluid–solid conjugate scalar transport, accommodating arbitrary jumps in scalar and its flux. It eliminates the need for explicitly body-fitted grids, simplifies mesh generation for complex geometries, and enables application of standard finite-volume tools. Relative errors typically remain below $3\%$ in challenging multiphase test cases.

A plausible implication is that this approach facilitates accurate and efficient simulation in cases with complex interface geometry and general jump conditions, which are prevalent in thermal and chemical engineering contexts. However, grid convergence is first-order rather than second-order in the vicinity of the interface, due to the diffuse treatment; this may limit its use in applications where higher-order interface capturing is essential (Liu et al., 15 Jan 2026).

7. Related Directions and Significance

The VPM methodology constitutes a significant development for simulating conjugate transport across multiphase domains with nontrivial boundary conditions. Its validation against reference body-fitted mesh methods, as well as its compact divergence-form source for interfacial jumps, position it as a promising approach for complex engineering scenarios without requiring elaborate mesh construction. The method’s utility for cases with both scalar and flux jumps, combined with robust error characteristics and absence of nonphysical interior contamination, underscores its practical value in computational studies of multiphase transport phenomena (Liu et al., 15 Jan 2026).