Repo - API! - Replace FourFitVars with immutable MatrixSplines split into ideal and kinetic

[jhalpern30]

Release note

• Audience: developers
• Numerical impact: none (harness @ bc15b72)
• Migration: ForceFreeStatesResult.ffit is now .mats. ForceFreeStates.make_matrix and make_kinetic_matrix are now build_matrix_splines and build_kinetic_matrix_splines. Matrix fields are renamed after what they hold: amats → A_spline, fmats_lower → F_spline_lower, paats → P_spline_adj, and so on.

Completes the struct cleanup of #139 for ForceFreeStates (FFS). The Fourier-fitted stability matrices now live in an immutable MatrixSplines container holding an IdealMatrices struct and, on kinetic runs only, a KineticMatrices struct. Whether a run is kinetic becomes a property of the data — mats.kinetic === nothing — rather than something each caller infers from control flags, and the two physics models' matrices can no longer be mixed by accident.

Regression report

Run at bc15b720 against the stack base refactor/forcefreestates-reorg (#400), which is the parent
of this branch. Comparing against develop would fold #400's file reorganization into the same
report and make any movement unattributable between the two PRs.

regress --cases diiid_n1,solovev_n1,diiid_n1_riccati,gal_resistive_pe,solovev_kinetic_calculated,solovev_kinetic_ntv \
        --refs origin/refactor/forcefreestates-reorg,local

Every tracked quantity is bit-identical — all diffs exactly 0.0e+00, all profile checksums
identical. Nothing moved.

Case	Covers	Result
diiid_n1	ideal forward + PerturbedEquilibrium (PE)	47 unchanged
solovev_n1	ideal forward, second equilibrium	21 unchanged
diiid_n1_riccati	Riccati propagators, Δ′ BVP, delta_coil	17 unchanged
solovev_kinetic_calculated	FKG kinetic path	14 unchanged
solovev_kinetic_ntv	PE + KineticForces (KF) torque quadrature	6 unchanged
gal_resistive_pe	Galerkin assembly/solve	8 N/A — see note below

solovev_kinetic_calculated is the load-bearing one: it is the only case that exercises the two
changes with real construction semantics rather than pure renaming — Kw_spline/Kt_spline are now
built inside _compute_fkg_matrices from the raw kw_flat/kt_flat, and itp_opts is re-derived
locally instead of carried on the struct.

Full report — diiid_n1 (47 quantities)

Regression Report: diiid_n1
=======================================================================================================================
Ref 1: origin/refactor/forcefreestates-reorg  @ 6e090741 (2026-08-18)
       env: julia 1.12.6, arm64-apple-darwin24.0.0, manifest ece8f5d1 (pinned), 6 threads/6 BLAS
Ref 2: local  @ local (2026-08-18)
       env: julia 1.12.6, arm64-apple-darwin24.0.0, manifest ece8f5d1 (pinned), 6 threads/6 BLAS
-----------------------------------------------------------------------------------------------------------------------
Quantity                                      origin/refactor/forcefreestates-reorg  local            Diff       Status
-----------------------------------------------------------------------------------------------------------------------
total energy Re(et[1])                        8.013943e-01                           8.013943e-01     0.0e+00    OK    
total energy Im(et[1])                        1.233500e-04                           1.233500e-04     0.0e+00    OK    
plasma energy Re(ep[1])                       -1.348114e+00                          -1.348114e+00    0.0e+00    OK    
vacuum energy Re(ev[1])                       2.149508e+00                           2.149508e+00     0.0e+00    OK    
vacuum matrix min eigenvalue                  1.873975e-01                           1.873975e-01     0.0e+00    OK    
plasma energy (all)                           [35 elem]                              [35 elem]        0.0e+00    OK    
vacuum energy (all)                           [35 elem]                              [35 elem]        0.0e+00    OK    
total energy (all)                            [35 elem]                              [35 elem]        0.0e+00    OK    
ODE steps (saved)                             2655                                   2655             0.0e+00    OK    
ODE steps (total)                             4738                                   4738             0.0e+00    OK    
q0                                            1.204212e+00                           1.204212e+00     0.0e+00    OK    
q95                                           4.781723e+00                           4.781723e+00     0.0e+00    OK    
beta_t                                        1.327024e-02                           1.327024e-02     0.0e+00    OK    
beta_n                                        1.372511e+00                           1.372511e+00     0.0e+00    OK    
internal inductance li1                       8.842230e-01                           8.842230e-01     0.0e+00    OK    
internal inductance li2                       7.080727e-01                           7.080727e-01     0.0e+00    OK    
internal inductance li3                       7.304309e-01                           7.304309e-01     0.0e+00    OK    
poloidal beta betap1                          6.680738e-01                           6.680738e-01     0.0e+00    OK    
poloidal beta betap2                          5.349836e-01                           5.349836e-01     0.0e+00    OK    
poloidal beta betap3                          5.518763e-01                           5.518763e-01     0.0e+00    OK    
# singular surfaces                           5                                      5                0.0e+00    OK    
singular psi locations                        [5 elem]                               [5 elem]         0.0e+00    OK    
singular q values                             [5 elem]                               [5 elem]         0.0e+00    OK    
current beta betaj                            4.236479e-01                           4.236479e-01     0.0e+00    OK    
plasma volume                                 1.829472e+01                           1.829472e+01     0.0e+00    OK    
plasma current                                1.152130e+00                           1.152130e+00     0.0e+00    OK    
mpert                                         35                                     35               0.0e+00    OK    
npert                                         1                                      1                0.0e+00    OK    
toroidal field bt0                            2.006573e+00                           2.006573e+00     0.0e+00    OK    
wall field bwall                              3.880145e-01                           3.880145e-01     0.0e+00    OK    
aspect ratio                                  2.845746e+00                           2.845746e+00     0.0e+00    OK    
elongation kappa                              1.708322e+00                           1.708322e+00     0.0e+00    OK    
q profile (checksum)                          ed7c21fd61df...                        ed7c21fd61df...  identical  OK    
pressure profile (checksum)                   e15550827bf1...                        e15550827bf1...  identical  OK    
Mercier D_I profile (checksum)                5a6fcb1c3a97...                        5a6fcb1c3a97...  identical  OK    
resistive interchange D_R profile (checksum)  6284a4c9a75a...                        6284a4c9a75a...  identical  OK    
ballooning Delta' profile (checksum)          44bf968c25d5...                        44bf968c25d5...  identical  OK    
island half-widths                            [5 elem]                               [5 elem]         0.0e+00    OK    
Chirikov parameter                            [5 elem]                               [5 elem]         0.0e+00    OK    
||resonant area-weighted field||              5.207739e-04                           5.207739e-04     0.0e+00    OK    
PE plasma energy                              3.422586e+00                           3.422586e+00     0.0e+00    OK    
PE vacuum energy                              3.174509e+00                           3.174509e+00     0.0e+00    OK    
PE surface energy                             5.841099e+00                           5.841099e+00     0.0e+00    OK    
PE toroidal torque                            -5.062793e-02                          -5.062793e-02    0.0e+00    OK    
NTV torque FGAR [N·m]                         5.587060e-01                           5.587060e-01     0.0e+00    OK    
NTV kinetic energy dW FGAR [J]                6.903492e-02                           6.903492e-02     0.0e+00    OK    
Runtime (s)                                   161.7s                                 162.4s                      --    
resonant area-weighted field b^r              [5 elem]                               [5 elem]         0.0e+00    OK    
=======================================================================================================================
Summary: 47 unchanged

Full report — solovev_kinetic_calculated (14 quantities)

Regression Report: solovev_kinetic_calculated
==============================================================================================
Ref 1: origin/refactor/forcefreestates-reorg  @ 6e090741 (2026-08-18)
       env: julia 1.12.6, arm64-apple-darwin24.0.0, manifest ece8f5d1 (pinned), 6 threads/6 BLAS
Ref 2: local  @ local (2026-08-18)
       env: julia 1.12.6, arm64-apple-darwin24.0.0, manifest ece8f5d1 (pinned), 6 threads/6 BLAS
----------------------------------------------------------------------------------------------
Quantity                 origin/refactor/forcefreestates-reorg  local          Diff     Status
----------------------------------------------------------------------------------------------
vacuum energy Re(ev[1])  1.037963e+01                           1.037963e+01   0.0e+00  OK    
q0                       1.900003e+00                           1.900003e+00   0.0e+00  OK    
total energy Re(et[1])   1.861536e+00                           1.861536e+00   0.0e+00  OK    
total energy Im(et[1])   -1.418013e+00                          -1.418013e+00  0.0e+00  OK    
plasma energy Re(ep[1])  -8.518098e+00                          -8.518098e+00  0.0e+00  OK    
ODE steps (saved)        600                                    600            0.0e+00  OK    
# singular surfaces      2                                      2              0.0e+00  OK    
singular psi locations   [2 elem]                               [2 elem]       0.0e+00  OK    
total energy (all)       [32 elem]                              [32 elem]      0.0e+00  OK    
ODE steps (total)        735                                    735            0.0e+00  OK    
mpert                    32                                     32             0.0e+00  OK    
q95                      3.147422e+00                           3.147422e+00   0.0e+00  OK    
Runtime (s)              132.9s                                 107.9s                  --    
singular q values        [2 elem]                               [2 elem]       0.0e+00  OK    
npert                    1                                      1              0.0e+00  OK    
==============================================================================================
Summary: 14 unchanged

Runtimes came in 5–20% lower on ref 2 across every case including the unchanged ones, which is
worktree cache ordering rather than anything this branch did. No performance claim is made.

Notes for reviewers

Reading order. The six commits are meant to be read in sequence, not squashed: freeze the flat
struct, split it into IdealMatrices/KineticMatrices, drop the helpers the split made redundant,
then three renaming passes. The first two carry the design; the rest are mechanical.

One rename was not mechanical. compute_node_xi_s! had a local named mats selecting the active
model. Renaming the ffit parameter to mats would have made that local shadow the parameter and
throw UndefVarError on its own right-hand side; the local is now active_mats, matching the name
FieldReconstruction.jl already used for the same selection.

A consistency guard was removed. CalculatedKineticMatrices.jl asserted
mats.numpert_total == ffs_intr.numpert_total. With numpert_total off the struct there is nothing
left to disagree — np now has a single source — but a mismatched mats would surface as a
dimension error inside the FKG loop rather than as a named assertion. Worth a look if you think that
guard was earning its keep.

gal_resistive_pe tracks nothing today, and that is pre-existing. All 8 quantities report N/A on
both refs, so it is not a regression from this branch. A confirming run of develop @ bb595659
against the stack base reports the same 8 N/A, so the case is dead on develop too — this is neither
#383 nor #400. The cause is in the example deck's own
header: "The PE stage currently warns and skips: it requires the free-boundary δW, which the
Galerkin formalism does not yet produce." The case tracks 8 PerturbedEquilibrium/SingularCoupling/*
datasets that the run never writes, so it costs ~4 minutes per harness invocation and guards nothing
until gal-side δW lands. Flagging it as harness debt, not as anything to fix here.

A second harness gap. No case exercises kinetic_source = "fixed" with kinetic_factor > 0 —
every example deck using the fixed source sets the factor to 0.0, and the only deck with a nonzero
factor uses "calculated". fixed_kinetic_matrices changed signature here (it now receives np
explicitly instead of reading it off the struct, including the np ÷ mpert multi-n tiling loop).
That path is covered by runtests_fullruns.jl at kinetic_factor = 1e-9 in both single-n and
multi-n, so it is verified — but by tests, not tracked by the harness.

A latent defect was found and deliberately not fixed here. compute_clebsch_displacements
(FieldReconstruction.jl) applies cholesky!(Hermitian(amat, :L)) to whichever A matrix the active
model supplies. In a kinetic run that is the kinetic A, which is not Hermitian. The defect predates
this PR — before the split, the same line read the same kinetic-overwritten matrix — so this branch
neither introduces nor fixes it, and the fix is a numerical change that does not belong in a PR
claiming inertness. It needs its own branch, its own harness run, and a regression case; the Fortran
reference (gpeq.f) turns out to skip the factorization entirely on the kinetic path, which is
probably the right shape for the fix.

[jhalpern30]

@matt-pharr I think this is ready to go. I will be adding one more PR to the stack that should be super small regarding the bug references in "A latent defect was found and deliberately not fixed here" above. I'm assigning you to it since I figured you would just merge whenever you're ready to continue your workflow, but lmk if you want me to fix anything. I am pretty sure we discussed everything I changed here, so should be a straightforward review

[jhalpern30]

I am going to also say this PR officially closes #139. All the discussion in that issue is now out of date anyway. The remaining related items would be extracting intr and odet into substructs, which we've also talked about