[Paper Review] Pick the Reconstruction with Smaller Jumps — Deng (2018) MUSCL-THINC-BVD

Choosing between two candidate reconstructions by minimizing total boundary variation

In 2018 the Xiao group at Tokyo Tech ran a measurement nobody likes to publish. They simulated a compressible two-phase flow for 1000 steps and measured the interface thickness afterward. The jump that started one cell wide had spread across eight cells. WENO did not help much either. Even schemes whose dissipation looks small per step accumulate into a smeared interface over long simulations. This post unpacks their answer to that problem — the MUSCL-THINC-BVD reconstruction. The idea is simple. For every cell, build two candidate reconstructions at once, then keep the one whose jumps at cell boundaries are smaller.

One-Page Summary#

Authors / Journal: Deng, Inaba, Xie, Shyue, Xiao. Journal of Computational Physics 371 (2018) 945–966.
Target problem: Material interfaces in the five-equation model for compressible two-phase flow get progressively smeared over time.
Proposal: At every cell, build a MUSCL reconstruction (for smooth regions) and a THINC reconstruction (for jumps). Compare the jump of each candidate at the cell boundaries and keep the smaller one. This is the BVD (Boundary Variation Diminishing) rule.
What's different: No post-processing anti-diffusion or artificial compression. The same BVD rule is applied to volume fraction and to every other conservative variable, so consistency between variables is automatic.

Two Demands That Conflict#

Compressible two-phase numerics asks for two things at once. Smooth regions need accuracy and low dissipation. Jumps (interfaces, shocks) need monotonicity and a thin profile. One function rarely does both jobs well.

MUSCL is a piecewise-linear reconstruction. It guarantees monotonicity but is only first-order accurate in derivative. Every time it crosses an interface it shaves a little off, and over many steps that erodes the jump. THINC fits a tanh inside each cell and keeps a step within one or two cells, but applied to a smooth region it manufactures fake staircases.

The paper's starting point is a decision: do not blend the two into a single function. For each cell, just choose one.

MUSCL — The Reliable Diffuser#

The baseline candidate is MUSCL with the minmod slope limiter.

\tilde{q}_i^{\,\text{MUSCL}}(x) = \bar{q}_i + \sigma_i (x - x_i), \quad \sigma_i = \text{minmod}(\bar{q}_i - \bar{q}_{i-1},\; \bar{q}_{i+1} - \bar{q}_i)

Here $\bar{q}_i$ is the cell average, $\sigma_i$ is the slope at the cell center, and minmod returns the smaller magnitude when the two differences have the same sign and zero otherwise.

The boundary values produced by this candidate:

q^{L,\text{MUSCL}}_{i+1/2} = \bar{q}_i + \tfrac{1}{2}\sigma_i, \quad q^{R,\text{MUSCL}}_{i-1/2} = \bar{q}_i - \tfrac{1}{2}\sigma_i

Safe everywhere. No oscillations, but the jump thickness slowly grows.

THINC — A Monotone Function That Mimics Jumps#

THINC (Tangent of Hyperbola for INterface Capturing) fits a hyperbolic tangent jump inside the cell.

\tilde{q}_i^{\,\text{THINC}}(x) = \bar{q}_{\min} + \tfrac{\bar{q}_{\max}}{2}\left[1 + \theta \tanh\!\left(\beta\!\left(\tfrac{x - x_{i-1/2}}{\Delta x} - \tilde{x}_i\right)\right)\right]

$\bar{q}_{\min}, \bar{q}_{\max}$ are the local min and amplitude taken from the neighbor cells, $\theta = \mathrm{sgn}(\bar{q}_{i+1}-\bar{q}_{i-1})$ encodes the jump direction, $\beta$ controls the jump thickness, and $\tilde{x}_i$ is the jump center solved so that the cell average is preserved.

A $\beta$ between 1.4 and 2.0 is stable; 1.6 is the standard value. Larger $\beta$ confines the jump to a single cell — great where there is one, awful where there is not.

BVD — Measuring the Variation Between Two Candidates#

Now the core idea. With two candidates available, ask which one produces smaller jumps when it meets the neighbor reconstructions at the cell boundaries. The Total Boundary Variation (TBV) for cell $i$ is defined as

\text{TBV}_i^{P} = |q^{L,\text{MUSCL}}_{i-1/2} - q^{R,P}_{i-1/2}| + |q^{L,P}_{i+1/2} - q^{R,\text{MUSCL}}_{i+1/2}|

where $P$ is the candidate (MUSCL or THINC). Both neighbors are held fixed at MUSCL, and only the middle cell switches. The selection rule:

\tilde{q}_i^{\,\text{BVD}}(x) = \begin{cases} \tilde{q}_i^{\,\text{THINC}}(x) & \text{if } \delta < C_i < 1 - \delta,\ (\bar{q}_{i+1}-\bar{q}_i)(\bar{q}_i-\bar{q}_{i-1}) > 0,\ \text{TBV}^{\text{THINC}}_i < \text{TBV}^{\text{MUSCL}}_i \\ \tilde{q}_i^{\,\text{MUSCL}}(x) & \text{otherwise} \end{cases}

$C_i = (\bar{q}_i - \bar{q}_{\min})/\bar{q}_{\max}$ measures where the jump sits, $\delta \approx 10^{-4}$ . The first two conditions check whether a jump is plausible at all; the third is the actual decision.

Consistency on the Five-Equation Model#

This is where the paper goes beyond a reconstruction comparison. The five-equation model carries the volume fraction $\alpha_1$ , phasic densities $\alpha_1\rho_1, \alpha_2\rho_2$ , momentum $\rho u$ , and total energy $E$ . The same BVD rule is applied to every one of those variables.

Earlier methods sharpened the volume fraction and then patched the other variables to suppress pressure oscillations at the interface. BVD instead forces all five variables to share the same decision per cell. If a cell is judged "interface", every variable uses THINC there; if "smooth", every variable uses MUSCL. Consistency between variables falls out automatically and no separate anti-diffusion step is needed.

The BVD Decision Tree in NumPy#

import numpy as np
 
def muscl_minmod_edges(q):
    """ Take (q[i-1], q[i], q[i+1]) and return MUSCL left/right edge values """
    qm, q0, qp = q
    a, b = q0 - qm, qp - q0
    if a * b <= 0.0:
        slope = 0.0
    else:
        slope = np.sign(a) * min(abs(a), abs(b))
    return q0 - 0.5 * slope, q0 + 0.5 * slope  # qL, qR
 
def thinc_jump_edges(q, beta=1.6, eps=1e-20):
    """ tanh-fit jump: left/right edge values of THINC """
    qm, q0, qp = q
    qmin = min(qm, qp)
    qmax = max(qm, qp) - qmin
    if qmax < 1e-12:
        return q0, q0
    theta = np.sign(qp - qm)
    C = (q0 - qmin + eps) / (qmax + eps)
    if C <= 1e-6 or C >= 1 - 1e-6:
        return q0, q0
    B = np.exp(theta * beta * (2.0 * C - 1.0))
    A = (B / np.cosh(beta) - 1.0) / np.tanh(beta)
    qR = qmin + 0.5 * qmax * (1.0 + theta * A)
    num = 1.0 + theta * A * np.tanh(beta) + theta * np.tanh(beta)
    den = 1.0 + A * np.tanh(beta)
    qL = qmin + 0.5 * qmax * (num / den)
    return qL, qR
 
def bvd_decide(q_window, beta=1.6, delta=1e-4):
    """
    q_window: length-5 array [q[i-2], q[i-1], q[i], q[i+1], q[i+2]]
    returns: (qL_i, qR_i, picked_thinc)
    """
    qm2, qm1, q0, qp1, qp2 = q_window
    monotone = (qp1 - q0) * (q0 - qm1) > 0.0
    qmin = min(qm1, qp1)
    qmax = max(qm1, qp1) - qmin
    C = 0.5 if qmax < 1e-12 else (q0 - qmin) / qmax
    eligible = monotone and (delta < C < 1.0 - delta)
 
    qL_M, qR_M = muscl_minmod_edges([qm1, q0, qp1])
    if not eligible:
        return qL_M, qR_M, False
 
    qL_T, qR_T = thinc_jump_edges([qm1, q0, qp1], beta=beta)
    # TBV with both neighbors held at MUSCL
    _, qR_left = muscl_minmod_edges([qm2, qm1, q0])
    qL_right, _ = muscl_minmod_edges([q0, qp1, qp2])
    tbv_M = abs(qR_left - qL_M) + abs(qR_M - qL_right)
    tbv_T = abs(qR_left - qL_T) + abs(qR_T - qL_right)
    if tbv_T < tbv_M:
        return qL_T, qR_T, True
    return qL_M, qR_M, False
 
# Test: a smooth Gaussian + a square jump
N = 80
x = (np.arange(N) + 0.5) / N
q = np.where((x >= 0.6) & (x <= 0.8), 1.0, 0.0)
q += 0.5 * np.exp(-((x - 0.25) / 0.04) ** 2)
 
picks = 0
for i in range(2, N - 2):
    _, _, used_thinc = bvd_decide(q[i-2:i+3])
    picks += int(used_thinc)
print(f"THINC pick rate: {picks}/{N-4}")
# Expected: a few cells around the jump only

In the smooth Gaussian region, minmod captures the curvature and BVD prefers MUSCL. Only two or three cells on each side of the square jump pick THINC. More than 95% of cells stick with MUSCL — BVD does not deform the underlying profile.

Who Wins Each Cell — See It Directly#

The simulation below lets you run a 1D scalar advection with three reconstructions in parallel. Each step, the strip below the canvas highlights the cells in which BVD picked THINC.

N = 80β = 1.60CFL = 0.40

MUSCL THINC BVDt = 0.00THINC picked: 0% of cells

Sweeping $\beta$ from 1.4 to 2.2 shows the THINC cells get sharper at the cost of saw-tooth artifacts on the smooth Gaussian. Somewhere around 1.6 sits the sweet spot — that is why the paper chose it. Watch the strip too: only a few cells around the square jump light up cyan, while the smooth bump never triggers THINC. That is BVD's self-localization in action.

How Far You Can Trust It#

Convergence rate is below second order. MUSCL caps the smooth-region accuracy. Swapping the candidate for WENO recovers fifth-order, but this paper sells simplicity.
The TBV definition assumes the neighbors are MUSCL. The case where two adjacent cells both want THINC is folded away. Single-sweep decision works well in the paper's tests but later work generalizes it.
In multi-D, $\beta$ depends on the face-normal direction. The interface normal is estimated by a Young's algorithm. At kinks and corners the normal estimate degrades, so the jump can lean.
Interaction with high-order time integration. The paper uses SSPRK3 and re-evaluates the BVD decision at every RK stage. If the decision flickers between stages, the interface can wobble — a known follow-up topic.

Three-Line Takeaway#

For every cell, build two candidates — MUSCL for smooth, THINC for jumps — and keep the one whose total boundary jump is smaller.
Apply the same BVD rule to every conservative variable. Variable consistency (no spurious pressure oscillations) follows automatically.
No anti-diffusion or artificial compression is needed, yet interfaces stay one cell thin. Simplicity is the strength.