[logan-nc]

Release note

• Audience: users
• Numerical impact: none in any tracked quantity (harness @ 618e48d): all 27 cases unchanged. No harness deck reaches the kinetic-crossing path (see the probes below); gal_resistive_pe is N/A by design.
• Migration: decks run unchanged. A kinetic run that reaches a kinetic singular surface right after a reduction now does the forced reduction there, as Fortran does. That reduction leaves the solution unchanged; it only records sing_flag and resets the reference norms, which can shift later renormalizations, so results move at most at round-off. A qlim below q_min, from qhigh or from truncate_at_rational_offset, now stops with an error that names the cause instead of a failed root find.

These are the ForceFreeStates items of #591, the Fortran GPEC 1.6.0 corrections:

• a singular surface always forces a Gaussian reduction, even right after another reduction, and that repeated reduction keeps the previous pivot order (port of Fortran PR 332);
• with truncate_at_rational_offset on and no rational surface left to truncate at, the run stops with a clear error instead of a failed root find (Julia-only check; Fortran PR 323's empty-surface fixes need no port, see the probes).

What changed

Commit	#591 item	Fortran 1.6.0	Change
93255e19a	A6	PR 332 + 1f8abc54	compute_solution_norms!: the early return odet.new && !at_ideal_crossing reset the reference norms and skipped the reduction that sing_flag should force. A kinetic crossing, which passes no resonant rows, reached right after a reduction never reduced. Now fix_required = sing_flag as in Fortran. When the call comes right after a reduction, unorm is first rebuilt from that reduction's odet.index, so the repeated reduction pivots in the same order. Before, it used raw norms.
b89935112	A6	—	Unit test: a singular call with new = true reduces (ifix advances, sing_flag recorded) in the previous pivot order [3, 1, 2], not the raw-norm order [1, 2, 3].
4be641761	A12	— (Julia-only)	sing_lim! errors when stepping qlim back to keep rational_offset_fraction lands below q_min. Fortran sing_lim has no such check: its Newton search on the extrapolated spline can return ψ_lim < 0 instead of stopping with "Can't find psilim".
29bdbaae3	A12	—	The q_min check now covers every path to the root find, including truncate_at_rational_offset = false with qhigh below q_min. The message no longer suggests raising qhigh, which does nothing when qmax or the layer-overlap cap is what binds. The test no longer depends on where q_min falls.
b181766f9	—	—	Formats the branch with the pinned hooks, on top of #356's sweep.

Probes (temporary logging, not committed):

• A6: none of the harness kinetic decks reach the bug. The Solovev kinetic examples have one kinetic surface with no reduction just before it, and the DIII-D kinetic deck has no kinetic singular surfaces at n = 1 or n = 1, 2. That is why the unit test builds the state directly.
• A12: a Solovev run with no n = 1 rational surface in the domain (q_min = 1.90, qhigh = 1.95) completes with truncate_at_rational_offset = false, so the empty surface list that Fortran PR 323 guarded is already safe in Julia: every next-surface search in chunk_el_integration_bounds tests ising <= msing (kmsing) before indexing. With truncation on, the run used to die in Roots.ConvergenceFailed.

Regression report

Harness on feynman: the formatting sweep 1883beeb0 (#356) vs this branch at 618e48dff
tj_n1                          25 unchanged
diiid_slayer_n1_toroidal       5 unchanged
gal_resistive_diiid_galinner   11 unchanged
ggj_ray_q500i                  4 unchanged
solovev_kinetic_nuzero         14 unchanged
diiid_error_field              10 unchanged
diiid_slayer_n1                22 unchanged
diiid_truncation_n1            12 unchanged
gal_resistive_pe               8 missing/N/A
tj_kinetic_calculated          20 unchanged
diiid_slayer_n1_coupled        7 unchanged
gal_resistive_diiid            11 unchanged
solovev_n1                     22 unchanged
ggj_reference                  4 unchanged
tj_kinetic_ntv                 7 unchanged
solovev_multi_n                15 unchanged
diiid_n1                       53 unchanged
diiid_slayer_n1_chi_perp_e     22 unchanged
diiid_n1_riccati               17 unchanged
lar_ggj_tearing_n1             3 unchanged
diiid_multi_n                  19 unchanged
diiid_slayer_n1_winding        22 unchanged
efit_fixedbdy_separatrix       5 unchanged
tj_kinetic_multiion            6 unchanged
diiid_kinetic_multi_n          8 unchanged
diiid_n1_riccati_precompiled   17 unchanged
lar_reversed_shear_n1          22 unchanged

The baseline is #356's sweep, which develop now contains (5503af4fa); that sweep is itself bit-identical to develop (#356's report). The head 2515fca8f only adds that develop merge on top of 618e48dff, with no src changes.

Notes for reviewers

• Base: develop, merged up to 5503af4fa.
• Why the Fortran damage cannot occur in Julia: in Fortran the bug either left no sing_flag record (so idcon.f never zeroed the resonant column) or, in the intermediate 87e6ad16, re-sorted on raw norms and replaced the wrong index(1) column, moving W_t by ~2%. Julia's ideal crossing always reduced and records zeroed_idx, and it picks the replaced columns by their resonant rows, not by index(1). Its kinetic crossing is the con_flag = true port, which uses neither. This PR is Fortran parity in control flow, consistent with the unchanged harness.
• Change on the ideal path: when an ideal crossing comes right after a reduction, its non-resonant columns now pivot in the previous reduction's order instead of raw-norm order, as Fortran does. The resonant pivots are unchanged. After the resonant pivots every remaining column has zero u₁ in the resonant rows, so this changes only the basis of the same span.
• One gap left from Fortran, on purpose: a forced reduction before any earlier reduction (ifix == 0) uses raw norms. Fortran would use the axis-initialization order, but Julia keeps no such slot.
• ForceFreeStates - BUGFIX! - Trigger the Gaussian reduction on the conditioning of U₁ and factor U₁ in the plasma response matrix #520 rewrites the reduction trigger in the same function and must rebase over this. It does not touch the early return, so the two changes compose.
• Tests: alfven, at the develop-merged tip before formatting: runtests_eulerlagrange 137/137 (plus 26, 38 and 2 in its other top-level sets), runtests_layer_overlap_fitzpatrick 81/81 (9 + 43 + 29), runtests_sing 90/90, runtests_solve_api 173/173.
• Reviewed by an independent model reviewer. It found no correctness defect; its two improvements are in 29bdbaae3.

🤖 Generated with Claude Code

https://claude.ai/code/session_01VdRsGGg4YRG6E9fgLGbpdF

[matt-pharr]

Thankfully, not all of these were even a problem in julia in the first place thanks to @jhalpern30's fantastic work dealing with the ODE integrator... However, I am confident that this is a faithful port of the mentioned fixes that went into fortran v1.6.0.