RasmahRasmah

Mesh Quality

println("Julia ", VERSION)
Julia 1.12.6

Overview

A tetrahedral mesh is only as good as its worst element. Slivers and inverted tetrahedra make stiffness matrices ill-conditioned and gradients meaningless, so Rasmah exposes per-element quality measures and mesh-health helpers. This tutorial measures them, sees how a distorted mesh degrades them, and shows how remesh_optimize uses the same threshold to decide when to re-freeze topology.

The per-element measures

For a tetrahedron with edge matrix B = [p₂−p₁, p₃−p₁, p₄−p₁], Rasmah provides six scalar measures:

  • signed_tet_volume — det(B)/6; negative means an inverted (flipped) element.
  • tet_scaled_jacobian — det(B) / (‖b₁‖ ‖b₂‖ ‖b₃‖), in (−1, 1]; 1 is a perfect right-corner tetrahedron, 0 a degenerate (flat) one.
  • tet_quality (mean-ratio) — 3 ∛(det²) / Σ‖bᵢ‖², in [0, 1]; 1 is a right-corner tetrahedron (∛(4/9) ≈ 0.79 regular), 0 degenerate.
  • tet_aspect_ratio — longest / shortest edge; 1 regular, →∞ as a sliver collapses.
  • tet_radius_ratio — 3 · inradius / circumradius, in [0, 1]; 1 regular, →0 a sliver.
  • tet_condition_number — σmax/σmin of B; a regular tetrahedron scores 2 (its edge vectors can’t all be orthogonal), →∞ degenerate.
using Rasmah

mesh = tetrahedralize_box(1.0, 1.0, 1.0, 2, 2, 2)
e = mesh.elements[:, 1]

println("signed volume     = ", signed_tet_volume(mesh.nodes, e))
println("scaled Jacobian   = ", tet_scaled_jacobian(mesh.nodes, e))
println("mean-ratio quality = ", tet_quality(mesh.nodes, e))
println("aspect ratio      = ", round(tet_aspect_ratio(mesh.nodes, e); digits = 4))
println("radius ratio      = ", round(tet_radius_ratio(mesh.nodes, e); digits = 4))
println("condition number  = ", round(tet_condition_number(mesh.nodes, e); digits = 4))
signed volume     = 0.020833333333333332
scaled Jacobian   = 0.40824829046386296
mean-ratio quality = 0.5
aspect ratio      = 1.7321
radius ratio      = 0.7174
condition number  = 4.0489

The three added measures answer “how bad is the worst one?” from different angles: a sliver has a huge aspect ratio, a radius ratio near 0, and a large condition number — while a regular tetrahedron scores 1 on aspect/radius ratio, 1/√2 ≈ 0.707 on the scaled Jacobian, and 2 on the condition number:

The tetrahedralize_box cube-split is uniform: every element has the same quality. Let’s confirm, and check the volume sum matches the box:

qs = [tet_quality(mesh.nodes, mesh.elements[:, j]) for j in axes(mesh.elements, 2)]
println("quality range: ", round.(extrema(qs); digits=6))
println("sum of signed volumes = ", volume(mesh), "  (box = 1.0)")
quality range: (0.5, 0.5)
sum of signed volumes = 1.0000000000000007 m^3  (box = 1.0)

Mesh-health helpers

The whole-mesh health is summarized by min_quality and a quadratic quality_penalty that accumulates (1 − q)² over elements. needs_remesh flags a mesh whose worst element falls below a threshold:

println("min quality     = ", min_quality(mesh))
println("quality penalty = ", round(quality_penalty(mesh); digits=6))
println("needs_remesh (threshold 0.2) = ", needs_remesh(mesh; threshold=0.2))
min quality     = 0.40824829046386296
quality penalty = 16.808164
needs_remesh (threshold 0.2) = false

quality_report summarizes the whole distribution of the scaled Jacobian in one call (min, a few percentiles, max, and mean):

report = quality_report(mesh)
println(report)
(min = 0.40824829046386296, p5 = 0.40824829046386296, p50 = 0.40824829046386296, p95 = 0.40824829046386296, max = 0.40824829046386296, mean = 0.40824829046386296)

A uniform box is healthy. The threshold matters when a mesh gets distorted — exactly what happens to a frozen marching-tetrahedra mesh as its parameters move.

Watching quality degrade under deformation

Recall the frozen-topology idea from the shape-optimization tutorial: topology is fixed, and only vertex positions follow the parameters. As a hole moves, the tets stretch. Let’s track the worst scaled Jacobian as a bracket’s hole grows:

sdf(r) = to_sdf(difference(box(4.0, 2.0, 0.5), cylinder(r, 1.0)))
lo, hi = bounds(sdf(0.4))
margin = 0.05 .* (hi .- lo)

fm = freeze_marching_tetrahedra(sdf(0.4), lo .- margin, hi .+ margin, 8)
meshat(r) = compact(deform(fm, sdf(r)))

for r in (0.4, 0.6, 0.8, 1.0, 1.2)
    m = meshat(r)
    println(
        "r = ", r,
        "   min scaled-Jacobian = ", round(minimum(
            tet_scaled_jacobian(m.nodes, m.elements[:, j]) for j in axes(m.elements, 2)
        ); digits=4),
        "   needs_remesh = ", needs_remesh(m; threshold=0.2),
    )
end
r = 0.4   min scaled-Jacobian = 0.0032   needs_remesh = true
r = 0.6   min scaled-Jacobian = -0.9923   needs_remesh = true
r = 0.8   min scaled-Jacobian = -0.9923   needs_remesh = true
r = 1.0   min scaled-Jacobian = -0.9927   needs_remesh = true
r = 1.2   min scaled-Jacobian = -0.9923   needs_remesh = true

As the hole eats into the plate, elements stretch, the worst scaled Jacobian drops, and needs_remesh flips from false to true — the signal remesh_optimize listens to before re-freezing a fresh topology.

Next steps