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.
Summary
_is_gap_spuriousconverts a normalized root to "kHz" by dividing bytauk, butQ/taukis an angular frequency in rad/s. The resulting value is compared againstgap_kHz_threshold, so the effective threshold is off by a factor of 2π (~6.28).src/Tearing/Dispersion/GrowthRateExtraction.jl:441:and line 444 for the next root.
Why
Q/taukis rad/s, not HzTwo independent sources agree:
Fortran GPEC,
slayer/params.f:47-50:Qconvis in seconds andQis dimensionless, soomega = Q/Qconvcarries whateveromegawas — andomegahere 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) definesQ = S^(1/3)·ω_E·τ_Hwith the perturbation taken asψ ∝ e^(iky − iωt), makingωangular by construction.So
imag(Q)/taukis in rad/s; multiplying by1e-3yields krad/s, not kHz. To get kHz the conversion needs a/(2π).This also means the
omega_Hz/gamma_Hzdataset names inTearing/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 HDF5unitsattribute instead.Impact
Contained but real.
_is_gap_spuriousonly gates the:gapwarning flag on an extracted root — it does not change anyQ_root,omega_Hz, orgamma_Hzvalue. 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_thresholdwas 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:
or rename the parameter to
gap_krad_per_s_thresholdand 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.