[StuartBenjamin]

Release note

• Audience: users
• Numerical impact: for n > 1, Galerkin Δ′ next to a singular surface moves by up to ~3 % (to the Riccati value where Galerkin converges). n = 1 is unchanged. (harness @ fc33044)
• Migration: none.

gal_make_grid! sizes the resonant and extension elements next to a singular surface in units of the singular factor m − nq ≈ −n q1 (ψ − ψ_s). It divides by |n q1| above a surface (GalerkinGrid.jl:101) but by |q1| below the next one (:133, marked "faithful to gal.f"). For n > 1 the elements below a surface are therefore n times wider than those above it. This PR uses n q1 on both sides, as PrincetonUniversity/GPEC#326 does for RDCON's gal.f.

To check that the symmetric choice is the right one, the size ratio of the elements below and above each surface was scanned on two DIII-D reconstructions at n = 2, with all element sizes (gal_dx0, gal_dx1, gal_dx2) shrunk together at gal_nx 512 (run on a branch whose GalerkinGrid.jl is identical to develop's, with the below-side factor set at run time). Innermost surface (q = 1.5); Riccati on the same input is an independent method that does not use the Galerkin grid:

ratio below/above	174956 at 1490 ms, dx × 1 / 0.5 / 0.25 / 0.125	201883 at 1628 ms, dx × 1 / 0.5 / 0.25 / 0.125
2 (develop)	−0.7622 / −0.7534 / −0.7535 / −0.7536	−0.6087 / −0.6076 / −0.6054 / −0.6009
√2	−0.7538 / −0.7481 / −0.7482 / −0.7483	−0.6038 / −0.6024 / −0.5995 / −0.5937
1 (this PR)	−0.7384 / −0.7419 / −0.7419 / −0.7420	−0.5977 / −0.5959 / −0.5924 / −0.5853
1/√2	−0.7194 / −0.7356 / −0.7356 / −0.7357	−0.5917 / −0.5896 / −0.5855 / −0.5771
1/2	−0.7036 / −0.7303 / −0.7303 / −0.7303	−0.5869 / −0.5846 / −0.5800 / −0.5708
Riccati (singfac_min 1e-4 / 3e-5 / 1e-5)	−0.7418 / −0.7419 / −0.7421	−0.5993 / −0.5994 / −0.5996

On 174956, Δ′ converges as the elements shrink at every ratio, but to a different value for each ratio, and only ratio 1 matches Riccati (−0.7420 against −0.7418 to −0.7421). On the other surfaces the change is ≤ 0.2 %. For n = 1 the two expressions are the same, so nothing changes there.

The scan also shows a weakness of the Galerkin method that this PR doesn't fix: the innermost Δ′ depends on the size ratio of the elements on the two sides of a surface, and on 201883 it does not converge as the elements shrink, at any ratio.

Regression report

regress --refs 9fb32c9a6,fc330443f (develop vs this head after merging develop, all 19 cases; julia 1.12.6, 8 threads): every tracked quantity is unchanged except SLAYER γ_Hz in diiid_slayer_n1, which differs by 0.196 Hz (0.02 %) against a 0.1 Hz noise threshold. That case runs integrator = "riccati", so the Galerkin grid this PR changes is never reached; its layer inputs all match exactly and only the AMR root extraction moved (Q_root by 1.2e-5, ω by 0.014 Hz), which is run-to-run noise. gal_resistive_pe is N/A on both refs, as on develop. The Galerkin cases are n = 1, where the two expressions are the same; the n = 2 effect is in the tables above.

---

🤖 Generated with Claude Code

[matt-pharr]

Thank you for the fix Stuart! This one has been on my back burner for a long time. This matches the mentioned Fortran PR.