InnerLayer.SLAYER - BUGFIX! - Normalize D by the electron diamagnetic share, not the temperature ratio

[d-burg]

Release note

• Audience: users
• Numerical impact: D_norm moves on every deck with T_i ≠ T_e, including the shipped diiid_slayer_n1 case (D_norm 20.4%, γ_Hz 9.8% on diiid_slayer_n1). The ion-sound-radius normalization had the electron and ion temperatures transposed. (harness @ 9b44f8c)
• Migration: re-run any SLAYER case whose profiles have T_i ≠ T_e — saved D_norm, Q_root, ω_Hz and γ_Hz from earlier runs are not comparable to new ones, so archived numbers and any figures built from them need regenerating. No API change. One behavioural change to watch for: a caller passing same-sign ω_*e/ω_*i (unphysical, but previously silently accepted) now raises an ArgumentError from slayer_parameters instead of producing a value.

D_norm used √(τ/(1+τ)) with τ = T_i/T_e. Fitzpatrick's definition is D = S^(1/3)·ι_e^(1/2)·d̂_β, and ι_e reduces to 1/(1+τ), not τ/(1+τ). The two agree only at T_i = T_e. GPEC already computes the correct ι_e and did not use it.

The defect

From TJ Documentation/Layer.tex:

eq.	statement
(228)	D = S^(1/3) · ι_e^(1/2) · d̂_β
(282)	ι_e = Q_e/(Q_e − Q_i)
(213)	ι = (T_e/T_i)·((1+η_e)/(1+η_i)), with ι_e = ι/(1+ι)

With η_e = η_i the third line gives ι = T_e/T_i = 1/τ, hence

ι_e = (1/τ)/(1 + 1/τ) = 1/(1+τ) = T_e/(T_e + T_i)

GPEC used τ/(1+τ) = T_i/(T_e + T_i) — the same expression with the two temperatures swapped. Equal at τ = 1, divergent in opposite directions otherwise.

Note that (282) is exactly the iota_e GPEC already computes and stores from ω_*e/ω_*i. The fix is to use it.

Why the computed iota_e and not 1/(1+τ)

1/(1+τ) is only the η_e = η_i limit. The full ι also carries (1+η_e)/(1+η_i), the ratio of the electron and ion pressure-gradient scale lengths — TJ's own text calls ι "the ratio of the electron to the ion pressure gradient at the rational surface". GPEC's iota_e = Q_e/(Q_e − Q_i), built from the actual ω_*e and ω_*i spline derivatives, captures that automatically. It is both TJ-faithful and strictly more general than either closed form. TJ itself hardcodes ι_e = 0.5 (Rational.cpp:198), so this is a case where GPEC can do better than the reference.

Impact on the shipped deck

examples/DIIID-like_ideal_example/TkMkr_D3Dlike_Hmode_kinetic.h5 has T_i ≠ T_e — 2534.73 vs 1952.59 eV at index 200, so τ ≈ 1.298:

	value
old factor √(τ/(1+τ))	0.752
new factor √ι_e (η-equal estimate 0.435)	0.660
D_norm ratio	≈ 0.88

So this is not latent — unlike the companion n-factor PR, it changes diiid_slayer_n1 today.

Guard

√ι_e requires ι_e > 0, which holds whenever ω_*e and ω_*i have opposite signs — the physical case, and the sign convention the slayer_parameters docstring already states. Same-sign inputs with |Q_i| > |Q_e| would make ι_e < 0; that now raises a clear ArgumentError rather than a bare DomainError, mirroring the existing degenerate-Q_e == Q_i guard immediately above it.

Validation

• runtests_slayer_params 49/49 — the one failure this produced was runtests_slayer_params.jl:88, the assertion that pins the definition of D_norm; it is updated to the corrected formula in this PR. That is a definition change, not a loosened tolerance.
• runtests_slayer_riccati 29/29 — the large-D/small-D branch selection assertions still hold. Worth knowing: that fixture passes same-sign ω_*e/ω_*i (1e4, 5e3), which is unphysical and gives iota_e = 2.0, so D_norm roughly doubles there. The branch inequalities survive because the fixture sits well inside its regime, but the fixture arguably wants fixing separately — flagging rather than quietly retuning it.
• Full suite julia -t 4 --project=. test/runtests.jl at a5dea9fed: 62 testsets, 2345 assertions, 0 failures
• Regression harness — below

Regression report

diiid_slayer_n1, develop (e09795119) vs this head (a5dea9fed), both sides --force, identical pinned environments (julia 1.11.6, manifest 7e5c34ad).

Regression Report: diiid_slayer_n1
===========================================================================================
Ref 1: e09795119  @ e09795119 (2026-08-27)
       env: julia 1.11.6, arm64-apple-darwin24.0.0, manifest 7e5c34ad (pinned), 8 threads/8 BLAS
Ref 2: 9b44f8cfc  @ 9b44f8cfc (2026-08-28)
       env: julia 1.11.6, arm64-apple-darwin24.0.0, manifest 7e5c34ad (pinned), 8 threads/8 BLAS
-------------------------------------------------------------------------------------------
Quantity                            e09795119  9b44f8cfc  Diff                Status       
-------------------------------------------------------------------------------------------
SLAYER surface indices              [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER poloidal m                   [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER toroidal n                   [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER minor radius rs              [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER r-based shear                [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER Lundquist S                  [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER D_norm                       [6 elem]   [6 elem]   5.587e-01 (20.44%)  ** CHANGED **
SLAYER P_perp                       [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER tauk                         [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER iota_e                       [6 elem]   [6 elem]   0.0e+00             OK           
SLAYER Q_root [2/1,3/1,4/1]         [3 elem]   [3 elem]   2.078e-03 (0.30%)   ** CHANGED **
SLAYER ω_Hz [2/1,3/1,4/1]           [3 elem]   [3 elem]   1.079e+00 (0.01%)   ** CHANGED **
SLAYER γ_Hz [2/1,3/1,4/1]           [3 elem]   [3 elem]   7.756e+01 (9.85%)   ** CHANGED **
SLAYER no_root flags [2/1,3/1,4/1]  [3 elem]   [3 elem]   0.0e+00             OK           
SLAYER enabled flag                 1          1          0.0e+00             OK           
Runtime (s)                         127.0s     133.1s                         --           
===========================================================================================
Summary: 4 changed, 11 unchanged

The report is self-verifying

The change is surgical, and the table proves it rather than asserting it:

• iota_e is 0.0e+00. It was already computed correctly on develop — this PR only starts using it. That is the single strongest line in the table: it separates "we changed which quantity D is built from" from "we changed iota_e".
• Every other layer input is bit-identical — rs, shear, Lundquist S, tauk, P_perp, the surface and mode-number lists. Nothing upstream of D moved.
• D_norm moves 20.4%, and γ_Hz follows at 9.8%. That is the physical consequence: the ion sound radius sets the drift-tearing layer response.

The γ shift is real, not solver noise. diiid_slayer_n1's γ line has a known reproducibility floor from the threaded root search — 0.076/0.122/0.145 Hz per surface (#418), which showed up as 0.02% in #431. Here it is 77.45 Hz, 9.84% — nearly three orders of magnitude above that floor.

On my pre-run estimate: I predicted ≈12% on D_norm from a single mid-radius point assuming η_e = η_i. The measured 20.4% is larger because τ varies across the six surfaces and the computed iota_e also carries the (1+η_e)/(1+η_i) factor that the closed form drops. The direction and order of magnitude held; the point estimate was not a prediction of the max.

Note

Harness re-run at 9b44f8cfc after the formula-citation comment restore touched src/ and staled the 8df7535e3 stamp (itself a re-run after an earlier comment-only commit). Each run reproduces the same numbers — the commits changed comments only — but the stamp must name a commit the harness actually ran at, so each src/ touch gets a real run rather than a re-pointed stamp.

Notes for reviewers

Opened as a draft. This one changes shipped numbers and wants a physics sign-off before it goes anywhere near ready.

Found while auditing whether the #431 shear sign was load-bearing elsewhere. Independent of the companion n-factor PR — different file, no textual conflict, and the missing n cancels in iota_e, so the iota_e this PR starts using is already correct on develop.

Cross-checks worth confirming against Park 2022 (docs/resources/2022-Park-Parametric dependencies...) and Burgess 2026, which I have not read for this: TJ is one source and the SLAYER lineage may have intended a different closure.

Suggested reviewer @matt-pharr (InnerLayer focus); assignee @d-burg. @logan-nc may want to see this too given it moves a harness case.

[jhalpern30]

Trusting your verification on the physics here, looks ok. Just clean up the LLM comments and its good to go

[jhalpern30]

No need for this comment - the error text is self explanatory

[jhalpern30]

This comment also feels kinda "Claude explaining code that it removed that quickly becomes out of date".

Not sure if this even needs a comment, or if it does, maybe just a reference to the relevant equation in Fitzpatrick's paper if desired