Skip to content

FFS - IMPROVEMENT - Cap the core knot density of the EL coefficient splines - #408

Open
logan-nc wants to merge 4 commits into
performance/consistent-surface-theta-parametrizationfrom
performance/decoupled-el-matrix-grid
Open

FFS - IMPROVEMENT - Cap the core knot density of the EL coefficient splines#408
logan-nc wants to merge 4 commits into
performance/consistent-surface-theta-parametrizationfrom
performance/decoupled-el-matrix-grid

Conversation

@logan-nc

@logan-nc logan-nc commented Aug 18, 2026

Copy link
Copy Markdown
Collaborator

Summary

The Euler-Lagrange coefficient splines (fmats/kmats/gmats and the primitives) inherited every knot of the equilibrium grid. A cubic spline's third-derivative jumps at knots scale as (node error)/Δψ³, so the equilibrium's near-axis packing (Δψ ~ 1e-6 at high mpsi) amplifies even tolerance-level node error (~1e-9, measured) into huge C² kinks — and the adaptive integrator's step size becomes slaved to the knot spacing. Measured directly: core jump magnitudes grow ~34× per mpsi doubling, and an (tol/J)^¼ step model reproduces the observed step-count ladder.

The coefficients are near-cylindrical in the core and do not need that packing. This PR builds their splines on a subset of the equilibrium grid with core density capped at Δψ ≥ 0.05·ψ below ψ = 0.1. Node values are unchanged — only interpolation density — which is why the physics moves at the 1e-7–1e-8 level.

Stacked on #398 (route (a), same investigation); diff shows only the cap once #398 lands.

Measured (DIII-D stripped decks, route-(a) base, mpsi 512/1024)

mψ=512 mψ=1024
accepted EL steps 2768 → 2188 4403 → 3034
step growth per doubling 1.59× → 1.39×
warm run 12.8 → 10.4 s (−19%)
et[1] relative change 3e-8 3e-8
Riccati Δ′ diagonal (all 5 surfaces) 4e-7

Not a tuned hack: the rule's basis and its generalization

Form: near the axis every component is a Frobenius power law in ψ; power laws are scale-free, so log-uniform sampling (Δψ ≥ c·ψ) resolves them at constant relative accuracy. Constant: cubic interpolation of ψ^p on a log-uniform grid errs by ~(p·c)⁴/384, so c = 0.05 resolves even the steepest spectrum component (p = m_max/2 = 11) to ~2e-4 — and the physics responds far below that because the steep components carry vanishing solution amplitude. Region: ends at min(0.1, innermost rational − RATIONAL_RES_RADIUS), and no knot inside a rational's resolution window is ever removed — guards for decks (e.g. higher n) whose rationals reach the core; verified no-ops on every current case.

Cross-equilibrium check (four equilibria, two construction paths, three grid families):

case knots EL steps et[1] change
DIII-D efit, traced, log-family m1024 1025→575 −31% 3e-8
tj_analytic_direct, traced, log-family m1024 1025→575 −35% 1.3e-9
LAR, inversion, packed auto m1024 1025→575 −8% exact to 8 digits
Solovev, traced, ldp m512/m1024 505→468 ~0 ~4.5e-7 abs (3e-5 rel via the ±10.4 → 0.0146 cancellation)

The LAR row is the anti-over-fit witness: clean geometry with no noise to exploit, and the cap is still harmless.

Two properties reviewers should know

  1. No-op on production decks. The two-pass auto grid's core spacing already satisfies the cap — verified: 287 → 287 knots on the DIII-D auto deck, and the full harness (diiid_n1, solovev_n1, diiid_n1_riccati, gal_resistive_diiid) reproduces EQUIL - IMPROVEMENT - Sample every flux surface at the same straight-fieldline angles #398's numbers exactly. The cap engages only on explicitly packed fine grids (large mpsi), which is precisely the issue Performance: Why do large equilibrium splines slow down the code? #376 scenario.
  2. A fixed cap beats general curvature selection — measured, not assumed. Removal-error knot selection (drop any knot whose interpolated value matches within δ for every element of all 12 matrices) was implemented and swept: it reaches the physics plateau only at δ=1e-7 where it keeps 504/513 knots and buys nothing, and breaks et[1]/Δ′ at any δ loose enough to matter. The core knots that are safe to remove fail the interpolation test (real coefficient curvature, but the solution components multiplying them vanish as ψ^|m|/2). Interpolation error is the wrong objective; the fixed mechanism-based cap is the right scope for the ideal path. Record: handoff/issue376/RESULTS.md §22–§23 (experiment branch).

Verification

⚠️ Requires third-party human review before merging — do not merge without an approving review.

🤖 Generated with Claude Code

logan-nc and others added 3 commits August 18, 2026 08:44
…plines

The Euler-Lagrange coefficient splines inherited every knot of the equilibrium grid.
A cubic spline's third-derivative jumps at knots scale as (node error)/dpsi^3, so the
equilibrium's near-axis packing (dpsi ~ 1e-6) amplifies even tolerance-level node
error into huge C2 kinks, and the adaptive integrator's step size becomes slaved to
the knot spacing -- measured directly: core jump magnitudes grow ~34x per mpsi
doubling while a (tol/J)^(1/4) step model reproduces the observed step-count ladder.

The coefficients are near-cylindrical in the core and do not need that packing. Build
their splines on a subset of the equilibrium grid with core density capped at
dpsi >= 0.05*psi below psi = 0.1; node values are unchanged, only knot density.

Measured on DIII-D stripped decks (route-a base, mpsi 512/1024): accepted EL steps
2768 -> 2188 and 4403 -> 3034, per-doubling growth 1.59x -> 1.39x, warm run -19% at
mpsi=512, with et[1] unchanged to 3e-8 relative and the Riccati BVP Delta-prime
diagonal unchanged to 4e-7 on all five surfaces.

Interim fixed-cap form; a matrix-curvature knot selection is planned to replace the
fixed rule, with this commit as the fallback.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
…es knots

On the production two-pass auto grid the cap is a no-op (its core spacing already
satisfies the density rule), so the unconditional message was noise.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
… protect rationals

Reframes the cap as log-uniform sampling of the Frobenius region: near the axis every
component is a power law in psi, and cubic interpolation of psi^p on a log-uniform
grid errs by ~(p*c)^4/384, so c = 0.05 resolves even the steepest spectrum component
(p = mmax/2) to ~2e-4 while the physics responds far below that.

Two generalization guards, motivated by cross-equilibrium testing: the capped region
now ends at the innermost rational surface when that sits inside psi = 0.1, and no
knot inside a rational's RATIONAL_RES_RADIUS window is ever removed -- preserving the
Delta'-stencil structure for decks (e.g. higher n) whose rationals reach the core.
Both guards are no-ops on every current case, verified: DIII-D m512 and Solovev m512
reproduce the previous cap's step counts and knot sets exactly.

Cross-equilibrium check of the rule itself: tj_analytic_direct m1024 (analytic,
traced) 1470 -> 959 steps at et[1] 1.3e-9; LAR m1024 (inversion path, clean geometry)
951 -> 875 at et[1] identical to 8 digits; Solovev ldp m512/m1024 ~unchanged steps at
~4.5e-7 absolute et[1] shift (a +-10.4 cancellation amplifies this to 3e-5 relative).

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
@logan-nc

Copy link
Copy Markdown
Collaborator Author

Review package (visual companion to this PR — decoupling diagram, knot-density and step-distribution before/afters, the Frobenius-scale basis for the cap, the cross-equilibrium generalization table, and the why-not-general-curvature-selection sweep):

📦 https://claude.ai/code/artifact/e1ed1ea9-1362-40bd-9f52-ccd386118e21

Self-contained page; complements rather than duplicates the description and diff. If the link does not resolve, ask @logan-nc to enable sharing on it.

…netic matrices

The expensive bounce-average kernel ran on every equilibrium knot, which both starves
kinetic decks of resolution (the shipped calculated example uses mpsi=16, carrying
percent-level errors on the B/C totals) and wastes kernel evaluations on fine grids.

Adds certified_kinetic_grid: seed with a coarse skeleton of the ideal coefficient-
spline knots (endpoints and rational windows always retained; the seed must not scale
with the equilibrium grid), then batched certify-or-refine rounds -- spline the kinetic
increments on the current knots, propose all uncertified midpoints, evaluate the whole
batch with the existing threaded kernel, and certify or insert. The certificate is the
max-element residual against tolerance * max|total| for every consumed family including
the adjoint combination kw3 - kt3, so tolerance is held on the ideal+kinetic totals the
solver consumes, never on increments in isolation. Spacing floors: the Frobenius cap in
the core, RATIONAL_RES_SPACING outside, and RATIONAL_RES_SPACING/4 inside rational
windows so narrow layers cannot hide behind the floor where they are expected.

Control: kinetic_grid_tol (default 0 = off; knob-off verified bit-identical). Measured
on the Solovev calculated case at tol=1e-3: m64 deck 65 -> 93 knots (certification adds
knots where the coarse grid under-resolves, matching the pre-registered expectation);
m256 deck 150 kernel evaluations vs 255 (-41%) with et[1] within 9.4e-5 of the full
grid; tol sweep 3e-3..1e-4 is plateau-flat. The kernel now accepts an arbitrary psi
list (psis kwarg), threading pattern unchanged.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant