Summary
s_compute_mthinc_normals in src/simulation/m_thinc.fpp computes interface normals using raw index-space differences:
nr_x = (alpha(j+1) - alpha(j-1)) * 0.5
nr_y = (alpha(k+1) - alpha(k-1)) * 0.5
nr_z = (alpha(l+1) - alpha(l-1)) * 0.5
No physical grid spacing (dx, dy, dz) is used.
Why this is correct on uniform grids
The MTHINC formulation works entirely in reference space — the unit cell ξ ∈ [−½, ½]^ndim. f_mthinc_volume_integral and f_mthinc_face_average integrate over reference coordinates with face positions fixed at ±½. The reference-space gradient of α is:
∂α/∂ξ ≈ (α(j+1) − α(j−1)) / 2
which is exactly what the code computes (the ×0.5 cancels on normalization). On a uniform grid, reference-space and physical-space normals coincide (up to a global scale factor), so the result is correct.
Why this is inaccurate on non-uniform grids
On a non-uniform Cartesian grid, the correct reference-space normal for cell (j,k,l) is proportional to:
( ∂α/∂x · Δx_j, ∂α/∂y · Δy_k, ∂α/∂z · Δz_l )
where ∂α/∂x ≈ (α(j+1) − α(j−1)) / (x_{j+1} − x_{j−1}). The code uses (α(j+1) − α(j−1)) / 2 instead, which is only correct when the stencil width equals the cell width (i.e., uniform spacing). On a stretched grid the normal direction can be wrong, proportional to the local grid non-uniformity.
Dividing by physical dx/dy/dz (as sometimes suggested) is not the correct fix — that gives the physical-space gradient, which is inconsistent with the reference-space face-average integrals. The correct fix is to properly account for the mapping between physical and reference coordinates when computing the gradient.
Impact
- Uniform grids: no impact.
- Mildly stretched grids: small accuracy loss in MTHINC normal direction.
- Strongly stretched grids (e.g., near-wall refinement): potentially significant error in interface orientation, degrading compression quality.
References
- B. Xie & F. Xiao, JCP 349 (2017) — THINC-QQ (physical-space formulation for unstructured grids)
- F. Xiao et al. (2011), Li & Xiao (2014) — structured-grid MTHINC reference-space formulations
Summary
s_compute_mthinc_normalsinsrc/simulation/m_thinc.fppcomputes interface normals using raw index-space differences:No physical grid spacing (dx, dy, dz) is used.
Why this is correct on uniform grids
The MTHINC formulation works entirely in reference space — the unit cell ξ ∈ [−½, ½]^ndim.
f_mthinc_volume_integralandf_mthinc_face_averageintegrate over reference coordinates with face positions fixed at ±½. The reference-space gradient of α is:which is exactly what the code computes (the ×0.5 cancels on normalization). On a uniform grid, reference-space and physical-space normals coincide (up to a global scale factor), so the result is correct.
Why this is inaccurate on non-uniform grids
On a non-uniform Cartesian grid, the correct reference-space normal for cell (j,k,l) is proportional to:
where
∂α/∂x ≈ (α(j+1) − α(j−1)) / (x_{j+1} − x_{j−1}). The code uses(α(j+1) − α(j−1)) / 2instead, which is only correct when the stencil width equals the cell width (i.e., uniform spacing). On a stretched grid the normal direction can be wrong, proportional to the local grid non-uniformity.Dividing by physical dx/dy/dz (as sometimes suggested) is not the correct fix — that gives the physical-space gradient, which is inconsistent with the reference-space face-average integrals. The correct fix is to properly account for the mapping between physical and reference coordinates when computing the gradient.
Impact
References