InnerLayer.SLAYER - BUGFIX - Use the shear magnitude in the layer timescales and widths

[d-burg]

Release note

• Audience: users
• Numerical impact: none on any existing result. abs(s) == s wherever the shear is positive, so every surface that produced a value before produces the same value; this only turns a crash into a result on reverse-shear surfaces. (harness @ 66cc5ed)
• Migration: none.

slayer_parameters used the signed r-based shear sval_r in five expressions that are magnitudes. On a reverse-shear surface this made the Alfvén time negative, hence a negative Lundquist number, and tauk = lu^(1/3)*tau_h threw a DomainError before any layer quantity was produced.

The five sites

| line | expression | why it needs |s| |
|---|---|---|
| chi_par_lmfp (×2) | (2 R0 vte)/(√π n s Wd) | a parallel diffusivity; must not change sign with shear direction |
| Wd | 1/√((rs/R0) s n) | an island width; sqrt of a negative argument throws |
| :lar critical-Δ | √((n s)/(R0 rs)) | same |
| tau_h | R0 √(μ₀ρ)/(n s bt) | a timescale; negative tau_h ⇒ negative lu ⇒ lu^(1/3) throws |

The layer timescales and widths depend on |dq/dr|, not its sign. The module's own diffusive-resistive width derivation already states the convention as (n|s|).

How it was found

Scanning DIII-D reconstructions for a separate piece of work. Shot 153072_3415 has q dipping to 1.855 off-axis (q0 = 2.364), which puts two q = 2 surfaces in the plasma with the inner one on the negative-shear branch. _find_rational_surfaces is deliberately reverse-shear-safe — it segments between q-extrema — so it finds both and hands them straight to the layer build, which then throws. Every SLAYER path is affected, not just the caller that surfaced it.

No shipped deck or harness case has reverse shear, which is why this has gone unnoticed.

Validation

• runtests_slayer_params 49/49, runtests_slayer_riccati 29/29, runtests_innerlayer 57/57 on this branch
• Reproduced the original failure and confirmed the fix scans 153072_3415 end to end (16 surfaces located, including both q = 2 branches)
• Regression harness, diiid_slayer_n1, develop d7e6e6fe vs this head 66cc5ed08, both sides --force:

SLAYER surface indices / m / n        0.0e+00   OK
SLAYER minor radius rs                0.0e+00   OK
SLAYER r-based shear                  0.0e+00   OK
SLAYER Lundquist S                    0.0e+00   OK
SLAYER D_norm                         0.0e+00   OK
SLAYER P_perp / tauk / iota_e         0.0e+00   OK
SLAYER Q_root [2/1,3/1,4/1]           5.1e-06   OK
SLAYER ω_Hz  [2/1,3/1,4/1]            1.9e-06   OK
SLAYER γ_Hz  [2/1,3/1,4/1]   1.276e-01 (0.02%)  ** CHANGED **
SLAYER no_root flags / enabled        0.0e+00   OK
Summary: 1 changed, 14 unchanged

The γ line is solver noise, not this change, and the report itself proves it: every layer input is bit-identical at 0.0e+00 — including sval_r, S and tauk, which are exactly the quantities abs() would have moved had it done anything on this deck. With identical inputs, a 0.13 Hz shift in the extracted root can only come from the root search. That magnitude sits inside the measured same-source reproducibility of the threaded search (0.076 / 0.122 / 0.145 Hz per surface, #418), which is above the case's current 1e-1 threshold — the subject of a separate threshold PR.

No shipped deck has a negative-shear rational surface, so zero change is the expected result and this is it.

Note

Opened as a draft. Reviewer suggestions: @logan-nc or @jhalpern30 — the substantive question is whether the shear sign is load-bearing anywhere else in the layer physics (the ω_* orientation being the obvious candidate), which I have not audited.

[github-actions]

This pull request is missing a reviewer.

If you are not ready to name them, mark this pull request as a draft.
docs/development/contributors.md suggests lead developers to ask.
Merging is not blocked here, but no pull request may be merged without human review.

[github-actions]

This pull request is missing a reviewer.

If you are not ready to name them, mark this pull request as a draft.
docs/development/contributors.md suggests lead developers to ask.
Merging is not blocked here, but no pull request may be merged without human review.

[d-burg]

@logan-nc @jhalpern30 small pre-existing bug, split out of #298 so it stands on its own.

slayer_parameters used the signed r-based shear in five expressions that are magnitudes, so a reverse-shear surface produced a negative Lundquist number and threw DomainError in lu^(1/3) before any layer quantity existed. Found on DIII-D 153072_3415, whose q dips to 1.855 off-axis.

abs(s) == s for positive shear, so nothing that works today changes value — it only turns a crash into a result. The harness bears that out: 14 of 15 quantities bit-identical at 0.0e+00, including sval_r, S and tauk, which are precisely what abs() would have moved had it done anything on that deck.

The question worth a physics opinion: I applied abs() on the grounds that these five are all magnitudes (a parallel diffusivity, an island width, an Alfvén time). If the shear sign is load-bearing somewhere else in the layer physics — the ω_* orientation being the obvious candidate — I have not audited that.

Draft for now; happy to mark ready once someone confirms the sign argument.