SLAYER gap_kHz_threshold compares a rad/s quantity against a kHz threshold (2π error) #372

Summary

_is_gap_spurious converts a normalized root to "kHz" by dividing by tauk, but Q/tauk is an angular frequency in rad/s. The resulting value is compared against gap_kHz_threshold, so the effective threshold is off by a factor of 2π (~6.28).

src/Tearing/Dispersion/GrowthRateExtraction.jl:441:

γ_idx = imag(sorted_roots[idx]) / tauk * 1e-3   # kHz

and line 444 for the next root.

Why Q/tauk is rad/s, not Hz

Two independent sources agree:

• Fortran GPEC, slayer/params.f:47-50:

  Qconv=lu**(1.0/3.0)*tau_h ! conversion to Qs based on Cole
  Q=Qconv*omega

  Qconv is in seconds and Q is dimensionless, so omega = Q/Qconv carries whatever omega was — and omega here is an angular frequency.

• Park 2022 PoP (docs/resources/2022-Park-Parametric dependencies of resonant layer responses across linear, two-fluid, drift-MHD regimes.pdf) defines Q = S^(1/3)·ω_E·τ_H with the perturbation taken as ψ ∝ e^(iky − iωt), making ω angular by construction.

So imag(Q)/tauk is in rad/s; multiplying by 1e-3 yields krad/s, not kHz. To get kHz the conversion needs a /(2π).

This also means the omega_Hz / gamma_Hz dataset names in Tearing/Roots/ are misnomers (they hold rad/s and 1/s respectively). Renaming them is a breaking output-schema change and is deliberately not proposed here — see the discussion in #364, which documents the true units in the HDF5 units attribute instead.

Impact

Contained but real. _is_gap_spurious only gates the :gap warning flag on an extracted root — it does not change any Q_root, omega_Hz, or gamma_Hz value. The consequence is that the γ-gap heuristic fires at ~6.28× the intended separation, so roots that should be flagged as suspiciously isolated are silently accepted. No regression-harness quantity moves.

Worth checking whether the default value of gap_kHz_threshold was tuned empirically against the current (rad/s-based) behaviour — if so, fixing the conversion without also retuning the default would change which roots get flagged.

Suggested fix

Either divide by 2π at the conversion:

γ_idx = imag(sorted_roots[idx]) / tauk / 2π * 1e-3   # kHz

or rename the parameter to gap_krad_per_s_threshold and document that it is angular. The first is preferable if the name is to stay meaningful to users; it needs the default re-examined in the same change.

Provenance

Found while reviewing #364 (self-describing HDF5 metadata), which audits the units of every tearing output dataset. Pre-existing on develop; not introduced by that PR.