ForceFreeStates - BUGFIX! - Default to the DCON fixed-axis initial condition by d-burg · Pull Request #440 · OpenFUSIONToolkit/GPEC

d-burg · 2026-09-02T19:38:12Z

Release note

Audience: users
Numerical impact: free-boundary energies move for every run that used the default axis condition — the fundamental et[1] by ≤ 5e-6 relative on the DIII-D-like cases and ≤ 4e-4 on the Solovev cases, ep[1]/ev[1] by ~2e-5, the full et/ep spectra by 3–9 % (the stiff non-resonant modes, which is the defect); Δ′ (Riccati BVP diagonal and raw matrix, Galerkin), ballooning, singular-surface and equilibrium quantities bit-identical; downstream SLAYER γ 0.02 %, PerturbedEquilibrium energies/fields ≤ 0.6 %, NTV torques 1.1–1.3 %; PE toroidal torque on diiid_n1 5.1e-2 → 1.2e-4 (harness @ c4276de) — 12 cases, report below
Migration: runs that must reproduce pre-existing free-axis results set fixed_axis = false under [ForceFreeStates]; everyone else changes nothing.

Free-boundary ideal δW from the ForceFreeStates (FFS) solver could report a spurious −10²…−10⁵ "instability" on well-behaved diverted equilibria. The cause is the default axis initial condition: the Euler-Lagrange integration started from a free-axis Frobenius state that does not reduce to the axis regularity condition at practical psilow, so one eigenvalue of the plasma response matrix W_p was wrong in every run and passed through a pole wherever it crossed zero. The default is now the DCON condition U₁ = 0, U₂ = I (ξ^ψ = 0 at the axis), which reproduces Fortran DCON's crit(ψ), W_p spectrum and et[1] on identical geqdsks. The Frobenius start remains available as fixed_axis = false.

The defect

On a 27-point n = 1, no-wall, diverted DIII-D-like ramp-up scan (TokaMaker geqdsks, set_psilim_via_dmlim = true, integrator = "riccati"), FFS reported et[1] ≈ −8e4 … −1.5e2 over a contiguous block q₀ ≈ 3.4–2.7 while every neighbouring point was stable at +1.5. Fortran DCON (GPEC/bin/dcon, identical mlow = −14, mhigh = 27, qlim = 7.2, psilim = 0.988) on the same files: stable, plasma matrix +0.25 … +1.6e4, no negative eigenvalue.

Signed smallest eigenvalue of W_p⁻¹ (DCON's crit) along ψ, flagged slice (left) and a clean neighbour (right). Black: Fortran DCON. Orange: julia_GPEC with the old default — displaced from the first stored step and never recovering; on the flagged slice it stays negative to the edge. Green: this PR — on top of DCON (median relative error 8e-5 over the 38 stored points).

Mechanism, established by measurement (read-only diagnostics; four independent investigations plus an adversarial pass):

stage	verdict
axis initial condition (`initialize_el_at_axis!`, `compute_axis_init`)	cause. Old default: U₂ = I, U₁ diagonal 0.11…0.23 for m = −14…−3, +6.1/−9.8/−13.7 for m = 3/4/5 at `psilow = 0.01` (docstring promised ~ψ_low^(\|m\|/2)). DCON: U₁ = 0.
integrator	exonerated — `"forward"` reproduces the same signed λ_min(S) as `"riccati"` to 5 significant figures (−5.979e-4)
tolerances, `ucrit`, `extended_precision_bvp`	exonerated — bit-identical
coefficient matrices (A, B, C, H, composite G)	exonerated — A, B agree with Fortran to 1e-6; G positive definite in both codes; 0.4 % apart in the pole direction
free-boundary assembly (`Free.jl` vs `free.f`)	exonerated — algebraically identical `wp = u₂ u₁⁻¹ / ψ₀²`
dmlim truncation, Δ′ BVP	not involved — Δ′ smooth across the pole

The single-variable intervention fixed_axis = true, nothing else changed:

flagged slice, psilim 0.98828	old default	this PR	Fortran DCON
signed λ_min(S), S = W_p⁻¹	−5.98e-4	+3.18e-3	+3.18e-3
stiffest eig(W_p) (m = −4…0 direction)	−8.55e4, one negative eigenvalue	+1.609e4, none	+1.609e4, none
`crit(psilim)`	−0.410	+2.179	+2.179
`et[1]`	−8.37e4	+1.505	+1.425 (see note)
Δ′ diagonal	−9.6637, −5.2578, −23.2488, −36.7497	identical to 6 decimals	—

On the clean neighbour the old default's stiff eigenvalue was +2.76e5 against DCON's +1.61e4 (17×): the eigenvalue is wrong everywhere, and the "pole" is merely where it crosses zero.

Left: et[1] over the 27-point scan, old default (orange) vs this PR (green). Right: the stiffest W_p eigenvalue on the flagged and clean slices, both settings, against Fortran DCON.

Full-scan sweep with the fix: 27/27 points, zero spurious values, et[1] equal to the previously e_v-filtered physical value at every point to ≤ 0.013, Δ′(2/1) changed by ≤ 1e-4 relative.

What changes

src/ForceFreeStates/CoreTypes.jl — fixed_axis::Bool = true; new docstring entry explaining both settings and the measured discrepancy of the Frobenius start.
src/ForceFreeStates/EulerLagrange.jl — branch comments updated; compute_axis_init docstring gains a warning that it does not deliver its documented limit at practical psilow. No code path changed.
src/ForceFreeStates/Riccati/Driver.jl — docstring updated.
test/runtests_fixed_axis.jl (new, wired into runtests.jl) — two full free-boundary runs on a new 129×129 synthetic diverted ramp-up fixture (test/test_data/regression_rampup_fixed_axis/, 303 kB, analytic profiles). With the default: et[1] > 0, no negative W_p eigenvalue, stiff W_p eigenvalue within 2 % of Fortran DCON's +1.609e4 (this fixture: +1.6090e4). With fixed_axis = false: the same three assertions as @test_broken, so a repaired compute_axis_init reports an unexpected pass. Δ′ diagonal asserted independent of the axis condition (rtol 1e-3). ~2 min.
docs/src/assets/fixed_axis_*.png — the two figures above.

Not in this PR (separate issues)

After the fix a uniform factor ≈ 1.056 remains between julia_GPEC's ep/ev/et and Fortran DCON's (e.g. et[1] 1.505 vs 1.425); the W_p spectrum agrees to ~1e-4, so this sits in the eigenproblem normalisation (Free.jl generalised pencil with the power norm vs free.f standard problem + ξ†Jξ/V′ division).
compute_axis_init itself: repair to the documented asymptotic limit, or demote. Suspects are the diagonal 2×2 truncation of A₀ and the regular-branch selection rule.
Crossing-injected non-Hermiticity in S (1e-3 at ideal singular-surface crossings vs 1e-11…1e-7 in the core).

Verification

Flagged and clean slices on this branch: numbers in the tables above (GeneralizedPerturbedEquilibrium.main on the staged decks; scripts and logs available on request).
Regression harness: regress --cases diiid_n1,diiid_n1_riccati,diiid_slayer_n1,efit_fixedbdy_separatrix,gal_resistive_diiid,gal_resistive_pe,solovev_kinetic_calculated,solovev_kinetic_multiion,solovev_kinetic_ntv,solovev_kinetic_nuzero,solovev_multi_n,solovev_n1 --refs 987e63f3b,c4276de63 (julia 1.11.6, arm64, manifest 7e5c34ad pinned, 8 threads / 2 BLAS). Summary per case:

case	summary	fundamentals moved	largest relative changes
`diiid_n1`	16 changed, 31 unchanged	total Re(et[1]) 8.013710e-01→8.013747e-01 ; total Im(et[1]) 4.144103e-05→8.242217e-05 ; plasma Re(ep[1]) -1.348365e+00→-1.348331e+00 ; vacuum Re(ev[1]) 2.149736e+00→2.149705e+00	PE toroidal torque 99.77%; ODE steps (saved) 4.27%; ODE steps (total) 3.61%
`diiid_n1_riccati`	6 changed, 11 unchanged	total Re(et[1]) 8.038607e-01→8.038645e-01 ; plasma Re(ep[1]) -1.344628e+00→-1.344593e+00 ; vacuum Re(ev[1]) 2.148488e+00→2.148458e+00	total energy (all) 2.59%; ODE steps (saved) 1.96%; ODE steps (total) 0.06%
`diiid_slayer_n1`	1 changed, 14 unchanged	—	SLAYER γ_Hz [2/1,3/1,4/1] 0.02%
`efit_fixedbdy_separatrix`	5 unchanged	—	—
`gal_resistive_diiid`	10 unchanged	—	—
`gal_resistive_pe`	8 missing/N/A	—	—
`solovev_kinetic_calculated`	7 changed, 7 unchanged	vacuum Re(ev[1]) 1.040320e+01→1.040466e+01 ; total Re(et[1]) 1.940488e+00→1.941020e+00 ; total Im(et[1]) -1.464711e+00→-1.462700e+00 ; plasma Re(ep[1]) -8.462708e+00→-8.463640e+00	total energy (all) 2.91%; ODE steps (total) 0.68%; ODE steps (saved) 0.17%
`solovev_kinetic_multiion`	5 changed, 1 unchanged	total Re(et[1]) 5.748525e-01→5.750592e-01	NTV total (D+T+imp+e) [N·m] 1.27%; NTV Deuterium [N·m] 1.13%; NTV Tritium [N·m] 1.13%
`solovev_kinetic_ntv`	3 changed, 3 unchanged	root-area-weighted total Re(et[1]) 5.748525e-01→5.750592e-01	NTV ψ quadrature evaluations 19.05%; NTV torque fgar [Re, Im] 1.13%
`solovev_kinetic_nuzero`	7 changed, 7 unchanged	vacuum Re(ev[1]) 1.068141e+01→1.068359e+01 ; total Re(et[1]) 2.385073e+00→2.385267e+00 ; total Im(et[1]) -1.984713e+00→-1.980744e+00 ; plasma Re(ep[1]) -8.296337e+00→-8.298325e+00	total energy (all) 2.91%; ODE steps (total) 0.82%; ODE steps (saved) 0.33%
`solovev_multi_n`	7 changed, 8 unchanged	total Re(et[1]) -1.906004e+00→-1.905962e+00 ; total Im(et[1]) -3.696819e-05→-3.696578e-05	plasma energy (all) 7.56%; ODE steps (saved) 3.03%; total energy (all) 2.39%
`solovev_n1`	10 changed, 12 unchanged	total Re(et[1]) 5.748525e-01→5.750592e-01 ; total Im(et[1]) -1.054421e-03→-1.422565e-05 ; plasma Re(ep[1]) -9.821799e+00→-9.821498e+00 ; vacuum Re(ev[1]) 1.039665e+01→1.039656e+01	plasma energy (all) 9.33%; total energy (all) 2.87%; ODE steps (saved) 2.28%

Δ′ is bit-identical wherever it is computed (diiid_n1_riccati: BVP diagonal and raw matrix, diff 0.0; gal_resistive_diiid: PEST3 Δ diagonal, ‖Δ′‖, coil block, inner-layer Δ all 0.0). The imaginary parts of et[1], which should vanish, move toward zero on solovev_n1 (−1.05e-3 → −1.4e-5). One line reviewers should look at: PE toroidal torque on diiid_n1 drops from 5.1e-2 to 1.2e-4 — the ideal perturbed-equilibrium response should carry no net torque, so this reads as the old axis condition having injected a spurious one, but I have not traced it.

All changed quantities (62 lines)

diiid_n1 — DIII-D-like equilibrium, n=1, ideal + perturbed equilibrium
Summary: 16 changed, 31 unchanged

total energy Re(et[1]) | 8.013710e-01 | 8.013747e-01 | 3.726e-06 (0.00%)
total energy Im(et[1]) | 4.144103e-05 | 8.242217e-05 | 4.098e-05 (98.89%)
plasma energy Re(ep[1]) | -1.348365e+00 | -1.348331e+00 | 3.389e-05 (0.00%)
vacuum energy Re(ev[1]) | 2.149736e+00 | 2.149705e+00 | 3.017e-05 (0.00%)
plasma energy (all) | [35 elem] | [35 elem] | 3.476e+02 (2.95%)
vacuum energy (all) | [35 elem] | [35 elem] | 3.331e+00 (0.16%)
total energy (all) | [35 elem] | [35 elem] | 3.472e+02 (2.58%)
ODE steps (saved) | 2647 | 2534 | 1.130e+02 (4.27%)
ODE steps (total) | 4715 | 4545 | 1.700e+02 (3.61%)
island half-widths | [5 elem] | [5 elem] | 2.700e-05 (0.09%)
Chirikov parameter | [5 elem] | [5 elem] | 6.584e-04 (0.08%)
||resonant area-weighted field|| | 5.188788e-04 | 5.196903e-04 | 8.115e-07 (0.16%)
PE plasma energy | 3.422683e+00 | 3.420577e+00 | 2.107e-03 (0.06%)
PE surface energy | 5.826403e+00 | 5.828561e+00 | 2.159e-03 (0.04%)
PE toroidal torque | 5.086535e-02 | 1.155847e-04 | 5.075e-02 (99.77%)
resonant area-weighted field b^r | [5 elem] | [5 elem] | 2.419e-06 (0.60%)

diiid_n1_riccati — DIII-D-like equilibrium, n=1, Riccati integrator Δ' matrix
Summary: 6 changed, 11 unchanged

total energy Re(et[1]) | 8.038607e-01 | 8.038645e-01 | 3.820e-06 (0.00%)
plasma energy Re(ep[1]) | -1.344628e+00 | -1.344593e+00 | 3.416e-05 (0.00%)
vacuum energy Re(ev[1]) | 2.148488e+00 | 2.148458e+00 | 3.034e-05 (0.00%)
total energy (all) | [35 elem] | [35 elem] | 3.472e+02 (2.59%)
ODE steps (saved) | 51 | 50 | 1.000e+00 (1.96%)
ODE steps (total) | 1667 | 1666 | 1.000e+00 (0.06%)

diiid_slayer_n1 — DIII-D-like H-mode equilibrium, n=1, SLAYER tearing-mode analysis (uncoupled per-surface, AMR, validity-gated)
Summary: 1 changed, 14 unchanged

SLAYER γ_Hz [2/1,3/1,4/1] | [3 elem] | [3 elem] | 1.276e-01 (0.02%)

efit_fixedbdy_separatrix — Fixed-boundary EFIT (box hugs LCFS), Newton separatrix-find regression
Summary: 5 unchanged

gal_resistive_diiid — DIII-D-like, n=1, RDCON outer-region Galerkin Δ′ with rpec coil columns (delta_coil)
Summary: 10 unchanged

gal_resistive_pe — DIII-D-like, n=1, DRIVEN/RPEC: gal-matched resistive ξ → PerturbedEquilibrium (coil-driven singular coupling)
Summary: 8 missing/N/A

solovev_kinetic_calculated — Solovev, n=1, kinetic_source=calculated, kinetic_factor=1.0, 50/50 D-T multi-ion kinetic matrices, mode-converged band (delta_m=8, mpert=32); pressure-consistent kinetic.dat (n·T tied to the Solovev P₀≈4.27e4 Pa).
Summary: 7 changed, 7 unchanged

vacuum energy Re(ev[1]) | 1.040320e+01 | 1.040466e+01 | 1.464e-03 (0.01%)
total energy Re(et[1]) | 1.940488e+00 | 1.941020e+00 | 5.316e-04 (0.03%)
total energy Im(et[1]) | -1.464711e+00 | -1.462700e+00 | 2.010e-03 (0.14%)
plasma energy Re(ep[1]) | -8.462708e+00 | -8.463640e+00 | 9.323e-04 (0.01%)
ODE steps (saved) | 598 | 597 | 1.000e+00 (0.17%)
total energy (all) | [32 elem] | [32 elem] | 4.754e+02 (2.91%)
ODE steps (total) | 731 | 736 | 5.000e+00 (0.68%)

solovev_kinetic_multiion — Solovev n=1, ideal + PE + multi-ion NTV (D + T + electron summed)
Summary: 5 changed, 1 unchanged

NTV total (D+T+imp+e) [N·m] | 1.769531e-04 | 1.747001e-04 | 2.253e-06 (1.27%)
NTV Deuterium [N·m] | 1.055493e-04 | 1.043514e-04 | 1.198e-06 (1.13%)
NTV Tritium [N·m] | 1.109007e-04 | 1.096526e-04 | 1.248e-06 (1.13%)
NTV electron [N·m] | -3.949690e-05 | -3.930393e-05 | 1.930e-07 (0.49%)
total energy Re(et[1]) | 5.748525e-01 | 5.750592e-01 | 2.067e-04 (0.04%)

solovev_kinetic_ntv — Solovev analytical equilibrium, n=1, ideal + PE + NTV torque quadrature
Summary: 3 changed, 3 unchanged

NTV torque fgar [Re, Im] | [1 elem] | [1 elem] | 2.396e-06 (1.13%)
NTV ψ quadrature evaluations | 630 | 750 | 1.200e+02 (19.05%)
root-area-weighted total energy Re(et[1]) | 5.748525e-01 | 5.750592e-01 | 2.067e-04 (0.04%)

solovev_kinetic_nuzero — Solovev, n=1, kinetic_source=calculated, nutype=zero (collisionless energy integral), 50/50 D-T multi-ion kinetic matrices, mode-converged band (delta_m=8, mpert=32). Reuses the calculated example deck via override; pressure-consistent kinetic.dat.
Summary: 7 changed, 7 unchanged

vacuum energy Re(ev[1]) | 1.068141e+01 | 1.068359e+01 | 2.182e-03 (0.02%)
total energy Re(et[1]) | 2.385073e+00 | 2.385267e+00 | 1.940e-04 (0.01%)
total energy Im(et[1]) | -1.984713e+00 | -1.980744e+00 | 3.968e-03 (0.20%)
plasma energy Re(ep[1]) | -8.296337e+00 | -8.298325e+00 | 1.988e-03 (0.02%)
ODE steps (saved) | 598 | 600 | 2.000e+00 (0.33%)
total energy (all) | [32 elem] | [32 elem] | 4.743e+02 (2.91%)
ODE steps (total) | 736 | 742 | 6.000e+00 (0.82%)

solovev_multi_n — Solovev analytical equilibrium, multi-n, ideal stability
Summary: 7 changed, 8 unchanged

total energy Re(et[1]) | -1.906004e+00 | -1.905962e+00 | 4.218e-05 (0.00%)
total energy Im(et[1]) | -3.696819e-05 | -3.696578e-05 | 2.411e-09 (0.01%)
plasma energy (all) | [70 elem] | [70 elem] | 4.718e+02 (7.56%)
vacuum energy (all) | [70 elem] | [70 elem] | 3.029e+00 (0.02%)
total energy (all) | [70 elem] | [70 elem] | 4.717e+02 (2.39%)
ODE steps (saved) | 33 | 34 | 1.000e+00 (3.03%)
ODE steps (total) | 1114 | 1115 | 1.000e+00 (0.09%)

solovev_n1 — Solovev analytical equilibrium, n=1, ideal stability
Summary: 10 changed, 12 unchanged

total energy Re(et[1]) | 5.748525e-01 | 5.750592e-01 | 2.067e-04 (0.04%)
total energy Im(et[1]) | -1.054421e-03 | -1.422565e-05 | 1.040e-03 (98.65%)
plasma energy Re(ep[1]) | -9.821799e+00 | -9.821498e+00 | 3.011e-04 (0.00%)
vacuum energy Re(ev[1]) | 1.039665e+01 | 1.039656e+01 | 9.442e-05 (0.00%)
plasma energy (all) | [32 elem] | [32 elem] | 4.717e+02 (9.33%)
vacuum energy (all) | [32 elem] | [32 elem] | 1.830e+00 (0.02%)
total energy (all) | [32 elem] | [32 elem] | 4.717e+02 (2.87%)
ODE steps (saved) | 394 | 385 | 9.000e+00 (2.28%)
ODE steps (total) | 626 | 617 | 9.000e+00 (1.44%)
ca_left (checksum) | 2e29a01cdae9... | 259883ce4dc4... | 0.000e+00 (0.00%)

Tests: julia --project=. test/runtests.jl runtests_fixed_axis.jl runtests_eulerlagrange.jl runtests_riccati.jl — Fixed-axis default: 6 pass, 3 broken (the recorded Frobenius defect), 1m56s; ODE Tests 91 pass; materialize_derivative_stores! 26 pass; Riccati Integration Tests 14 pass.

🤖 Generated with Claude Code

…ndition `fixed_axis` now defaults to `true`: the Euler-Lagrange integration starts from U1 = 0, U2 = I at psilow (xi^psi = 0 at the magnetic axis), the condition used by Fortran DCON (dcon/ode.f, ode_axis_init) and Glasser 2016 Sec. VI. The previous default (`false`) started from the free-axis Frobenius state of `compute_axis_init`. At practical psilow (1e-2 ... 1e-4) that state has U1 entries of order 0.1-10 with mixed signs instead of the documented psilow^(|m|/2) limit. Started from it, the signed critical eigenvalue of W_p^-1 is displaced from the first stored step onward and never recovers: one eigenvalue of the plasma response matrix W_p sits ~10x off its DCON value in every equilibrium tested, and where it crosses zero the free-boundary energy et[1] reports a spurious -1e2 ... -1e5 "instability". Verified against Fortran DCON on identical geqdsks (diverted DIII-D-like, n = 1): with the fixed-axis start crit(psi) agrees to a median relative error of 8e-5, the stiff W_p eigenvalue to ~1e-4, and et[1] returns to the physical fundamental. The Delta' BVP is unaffected (<= 1e-4 relative). Integrator (riccati/forward), tolerances, coefficient matrices and the free-boundary assembly were all exonerated by measurement. The Frobenius start is retained behind `fixed_axis = false` for comparison; its docstring now records the measured discrepancy. Two figures documenting the defect and the fix are added under docs/src/assets/. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>

…iverted fixture Adds test/runtests_fixed_axis.jl and a 129x129 synthetic DIII-D-like ramp-up fixture (q0 = 3.42, diverted, no wall, n = 1; analytic profiles, 303 kB) on which the Frobenius free-axis start produces et[1] ~ -1e5 and a negative plasma-matrix eigenvalue while the DCON fixed-axis start gives et[1] = +1.50 and a stiff W_p eigenvalue of +1.609e4, matching Fortran DCON v1.5.5 on the 257x257 parent. Asserts: the default is fixed_axis = true; with it et[1] > 0, W_p has no negative eigenvalue and its stiff eigenvalue is within 2 % of the DCON value; the Delta' diagonal is independent of the axis condition (rtol 1e-3). The defect of the Frobenius start is pinned with @test_broken so a repaired compute_axis_init reports an unexpected pass. Runtime ~2 min (two full free-boundary runs at mpsi 128). Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>

Copilot

Warning

Copilot couldn't run its full agentic review because it didn't start before the timeout. Make sure your repository has a runner available, or add a copilot-code-review.yml file specifying one with the runs-on attribute. See the docs for more details.

Pull request overview

This PR changes ForceFreeStates’ default Euler–Lagrange axis initial condition to the DCON fixed-axis condition (fixed_axis = true) to eliminate spurious free-boundary instabilities and align crit(ψ), W_p, and et behavior with Fortran DCON.

Changes:

Default ForceFreeStatesControl.fixed_axis to true and document both axis-init options and their measured discrepancy.
Update related docstrings/comments in Euler–Lagrange and Riccati driver code paths to reflect the new default.
Add a new regression test and synthetic fixture to lock in the corrected default behavior and preserve the known defect as @test_broken for fixed_axis = false.

Reviewed changes

Copilot reviewed 8 out of 11 changed files in this pull request and generated 4 comments.

Show a summary per file

File	Description
`src/ForceFreeStates/CoreTypes.jl`	Switches `fixed_axis` default to `true` and adds detailed docstring rationale/impact.
`src/ForceFreeStates/EulerLagrange.jl`	Updates axis-init documentation and clarifies fixed-axis vs Frobenius behavior.
`src/ForceFreeStates/Riccati/Driver.jl`	Updates documentation to match the new default and DCON alignment claims.
`test/runtests_fixed_axis.jl`	Adds regression coverage for the default fixed-axis condition + records Frobenius defect as broken tests.
`test/runtests.jl`	Wires the new regression test into the test suite.
`test/test_data/regression_rampup_fixed_axis/gpec.toml`	Adds a test deck used by the new fixed-axis regression test.
`test/test_data/regression_rampup_fixed_axis/forcing.dat`	Adds forcing input for the new synthetic regression fixture.
`test/test_data/README_test_data.md`	Documents the new regression fixture and expected behaviors.

💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

+    wp = Matrix(fb.wp)
+    wpe = eigvals(Hermitian((wp + wp') / 2))
+    dp = real.(diag(Matrix(r.ffs.delta_prime.matrix)))
+    return (et1 = minimum(real.(fb.et)), wp_min = minimum(wpe), wp_max = maximum(wpe),


+    dir = mktempdir()
+    for f in readdir(FIXED_AXIS_FIXTURE)
+        cp(joinpath(FIXED_AXIS_FIXTURE, f), joinpath(dir, f))
+    end
+    toml = joinpath(dir, "gpec.toml")
+    s = read(toml, String)
+    s = replace(s, r"^fixed_axis.*\n"m => "")
+    s = replace(s, "[ForceFreeStates]\n" => "[ForceFreeStates]\nfixed_axis = $(fixed_axis)\n")
+    write(toml, s)
+    r = GeneralizedPerturbedEquilibrium.main([dir])
+    fb = r.ffs.free_boundary
+    wp = Matrix(fb.wp)
+    wpe = eigvals(Hermitian((wp + wp') / 2))
+    dp = real.(diag(Matrix(r.ffs.delta_prime.matrix)))
+    return (et1 = minimum(real.(fb.et)), wp_min = minimum(wpe), wp_max = maximum(wpe),
+            n_neg = count(<(0), wpe), dprime = dp)


+# Deprecated + ignored upstream (warning asks for removal), deleted here:
+#   mer_flag, ode_flag, mat_flag, use_parallel, force_wv_symmetry.


+# ForceFreeStates BVP (dW) + Delta'.
+#
+# Wall: nowall (CONVENTIONS.md §6 — required for valid Delta').
+# use_parallel: true  (CONVENTIONS.md §7 — false silently invalidates Delta').


d-burg and others added 2 commits September 2, 2026 15:16

d-burg self-assigned this Sep 2, 2026

github-actions Bot added bugfix Something was wrong and now is not changed-results Results move or an interface breaks - read before upgrading labels Sep 2, 2026

d-burg marked this pull request as ready for review September 2, 2026 22:37

d-burg requested review from logan-nc and a lite review from Copilot September 2, 2026 22:37

Copilot started reviewing on behalf of d-burg September 2, 2026 22:49 View session

Copilot AI reviewed Sep 2, 2026

View reviewed changes

Navigation Menu

Uh oh!

ForceFreeStates - BUGFIX! - Default to the DCON fixed-axis initial condition - #440

ForceFreeStates - BUGFIX! - Default to the DCON fixed-axis initial condition#440
d-burg wants to merge 2 commits into
developfrom
bugfix/forcefreestates-fixed-axis-default

d-burg commented Sep 2, 2026 •

edited

Loading

Uh oh!

Copilot AI left a comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

		# Deprecated + ignored upstream (warning asks for removal), deleted here:
		# mer_flag, ode_flag, mat_flag, use_parallel, force_wv_symmetry.

Uh oh!

Conversation

d-burg commented Sep 2, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Release note

The defect

What changes

Not in this PR (separate issues)

Verification

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Pull request overview

Reviewed changes

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

d-burg commented Sep 2, 2026 •

edited

Loading