Equilibrium - BUGFIX - Support original and rho grids for inverse equilibria

[StuartBenjamin]

Release note

• Audience: users
• Numerical impact: none for runs that worked before (solovev_n1 unchanged). grid_type = "original" and "rho" errored before for inverse equilibria. An inverse run with psilow = 0 on ldp, pow1, uniform or rho, which gave NaN at the axis, now gets an extrapolated axis row (harness @ f0497d3)
• Migration: none

Inverse equilibria now accept grid_type = "original" (the input file's own ψ nodes) and "rho" (uniform in √ψ).

Regression report

================================================================
Case: solovev_n1 — Solovev analytical equilibrium, n=1, ideal stability
================================================================

Regression Report: solovev_n1
=================================================================================
Ref 1: bfc60610e  @ bfc60610e (2026-10-09)
       env: julia 1.11.7, x86_64-linux-gnu, manifest fdfdedfa (pinned), 2 threads/2 BLAS, jit build
Ref 2: f0497d3f7  @ f0497d3f7 (2026-10-13)
       env: julia 1.11.7, x86_64-linux-gnu, manifest fdfdedfa (pinned), 2 threads/2 BLAS, jit build
---------------------------------------------------------------------------------
Quantity                      bfc60610e        f0497d3f7        Diff       Status
---------------------------------------------------------------------------------
total energy Re(et[1])        3.429538e-01     3.429538e-01     0.0e+00    OK    
total energy Im(et[1])        -7.747105e-06    -7.747105e-06    0.0e+00    OK    
plasma energy Re(ep[1])       -1.005916e+01    -1.005916e+01    0.0e+00    OK    
vacuum energy Re(ev[1])       1.040211e+01     1.040211e+01     0.0e+00    OK    
vacuum matrix min eigenvalue  2.171292e+00     2.171292e+00     0.0e+00    OK    
plasma energy (all)           [32 elem]        [32 elem]        0.0e+00    OK    
vacuum energy (all)           [32 elem]        [32 elem]        0.0e+00    OK    
total energy (all)            [32 elem]        [32 elem]        0.0e+00    OK    
ODE steps (saved)             409              409              0.0e+00    OK    
ODE steps (total)             407              407              0.0e+00    OK    
q0                            1.900006e+00     1.900006e+00     0.0e+00    OK    
q95                           3.147422e+00     3.147422e+00     0.0e+00    OK    
beta_t                        4.621020e-02     4.621020e-02     0.0e+00    OK    
beta_n                        3.215791e+00     3.215791e+00     0.0e+00    OK    
# singular surfaces           2                2                0.0e+00    OK    
singular psi locations        [2 elem]         [2 elem]         0.0e+00    OK    
singular q values             [2 elem]         [2 elem]         0.0e+00    OK    
mpert                         32               32               0.0e+00    OK    
npert                         1                1                0.0e+00    OK    
q profile (checksum)          90ffc338a738...  90ffc338a738...  identical  OK    
pressure profile (checksum)   5040338aec07...  5040338aec07...  identical  OK    
ca_left (checksum)            6f1c35be5425...  6f1c35be5425...  identical  OK    
Runtime (s)                   208.9s           217.5s                      --    
=================================================================================
Summary: 22 unchanged

Problem

For inverse equilibria (CHEASE, LAR), grid_type = "original" set mpsi but never built the ψ grid, so the run stopped with an undefined-variable error. Fortran inverse.f also offers rho, which Julia GPEC rejected as "Unsupported grid_type".

Fix

• original uses the input file's own ψ nodes, as Fortran does, and ignores psilow and mpsi (psihigh still sets the integration limit, as in Fortran); no minimum node spacing is applied. (Fortran takes the 1D input's nodes; they differ from the 2D ones only for efit_by_inversion.)
• rho is uniform in √ψ from psilow to psihigh: ψᵢ = (√psilow + (√psihigh − √psilow)·i/mpsi)², i = 0…mpsi. Fortran's rho ignores psilow (it runs from psihigh/(mpsi+1)² to psihigh); here it honours psilow like the other grids.
• When the first ψ node is the magnetic axis (the original grid of a file that starts at the axis, or psilow = 0), the straight-field-line weights and the GSE flux terms are 0/0 there. That row is now extrapolated from the three innermost surfaces, as deta already was: θ-dependent fields in √ψ (∝ r) by a new helper extrap_axis_row! (also used in equilibrium_gse!), dV/dψ in ψ.
• The grid_type docstring lists the new options.

Tests

New testset "Inverse grid types (original, rho)" in test/runtests_equil.jl, on the CHEASE test file INP1_binary (which starts at the axis): the original grid equals the file's nodes, all splines are finite and the extrapolated axis |B| is within 1 % of F/R₀ for every θ; the rho grid runs from psilow to psihigh and its splines are finite. test/runtests_equil.jl: 301/301 pass; the full suite was not run. Harness: solovev_n1 only, 22 unchanged (no harness case uses an inverse equilibrium on these grids; diiid_n1 was not run).

Notes for reviewers

#511 (field-line q for i-file inverse equilibria) edits the same sq_fs lines in equilibrium_solver; whichever lands second needs a small rebase.

---

🤖 Generated with Claude Code

[StuartBenjamin]

Superseded by #523 (same grids, without the axis-extrapolation code).