Build - Bump FastInterpolations.jl to v0.4.20 for PeriodicBC fix #531

Activity

1. StuartBenjamin

   assigned

   logan-nc

   on Oct 15, 2026
2. StuartBenjamin

   added

   bugfixSomething was wrong and now is notSomething was wrong and now is not

   on Oct 15, 2026
3. mgyoo86 commented on Oct 15, 2026

   

   mgyoo86

   on Oct 15, 2026

   Collaborator

   More actions

   Hi @StuartBenjamin,

   Thanks again for the this fix on FastInterpolations side :)

   The followings are some additional info for Logan and reviewers:
   this is a fairly small bug. Even on strongly non-uniform grids, the added interpolation error is around 1e-6, so it shouldn't have much impact on GPEC in practice. I also confirmed that the full test suite passes locally with v0.4.20.

   A bit more detail:

   • Where it applies. Most PeriodicBC call sites use uniform θ grids (range / rzphi_ys). There h_{i-1} = h_i, so the fix changes nothing. The main exception is DirectEquilibrium.jl:534, which splines the field-line ODE step points onto the uniform θ grid. Those points are somewhat irregular, mostly because of a short final step.
   • DIII-D-like example: δW, Δ′ and the equilibrium scalars move by ~1e-5 relative with a fixed ψ-grid, and by ~1e-4 with the default auto grid. Most of the auto-grid shift comes from the grid now picking 320 surfaces instead of 598.
   • Solovev example: ~1e-4–1e-3 at the default etol = 1e-10. At the example's loose etol = 1e-7 (only ~20–30 ODE points per surface), a few near-marginal quantities such as et[1] and Δ′ shift by up to a few percent.
   • In all cases v0.4.20 is closer to a tight-tolerance reference, and the stability verdicts are unchanged.

   v0.4.20 also brings other improvements, including performance work, so I also recommend updating the Project.toml to use FastInterpolations >= v0.4.20.
4. logan-nc commented on Oct 16, 2026

   

   logan-nc

   on Oct 16, 2026

   Collaborator

   More actions

   @mgyoo86 are there any API changes we need to be aware of? Or is it a fairly simple mater of just updating the Project.toml?
5. mgyoo86 commented on Oct 16, 2026

   

   mgyoo86

   on Oct 16, 2026

   Collaborator

   More actions

   @logan-nc
   This is just a minor/patch version bump under the same major version (v0.4.x), so there is no breaking changes for API at all.
   Therefore, updating the Project.toml is sufficient.
6. logan-nc

   added the

   Bug

   issue type on Oct 16, 2026
7. logan-nc

   added this to the GPEC v2.0.0 milestone on Oct 16, 2026
8. matt-pharr commented on Oct 16, 2026

   

   matt-pharr

   on Oct 16, 2026

   Collaborator

   More actions

   @mgyoo86 @jhalpern30 — I ran the full regression harness for the v0.4.20 bump (#532) and the shifts are much larger than the ~1e-5–1e-4 estimated above, so I bisected them. Summary: the PeriodicBC fix is correct and is the only part of v0.4.20 that matters, but its effect on Solovev is large, and it lands the DIII-D-like example's auto ψ-grid on a grid where the forward integrator breaks down.

   Harness, develop @ 278057d, manifests differing only in FastInterpolations (0.4.18 vs 0.4.20): 165 tracked quantities changed across 15 of 18 cases. Headlines:

   case	quantity	0.4.18	0.4.20
   diiid_n1 (forward, auto grid)	Re(et[1])	0.8009	1.7814
   	Im(et[1])	7.9e-5	0.531
   	ODE steps (saved)	1441	822
   diiid_n1_riccati (same eq. + grid)	Re(et[1])	0.8035	0.8032
   solovev_n1	Re(et[1])	0.442	0.0475
   equilibrium scalars (q0, q95, β, li)			~1e-6 rel.

   Bisect (examples run directly, develop source):

   FastInterpolations	Solovev et[1]	DIII-D-like et[1]	DIII-D ODE steps
   0.4.18	0.44197	0.80087	1441
   0.4.19	identical to 0.4.18	identical	1441
   0.4.20	0.04754	1.7814 + 0.531i	822
   0.4.20 with the du[i-1] line reverted	0.44197	0.80087	1422
   0.4.18 with only the du[i-1] fix	0.04754	0.80062	822

   So the one-line super-diagonal fix accounts for the shift; everything else in 0.4.19/0.4.20 is roundoff-level.

   Solovev — the shift is real and v0.4.20 converges correctly. Fixed ldp grid, tightening the equilibrium etol:

   etol	0.4.18	0.4.20
   1e-7 (example default)	0.442	0.0475
   1e-10	0.01606	0.01089
   1e-12	0.01603	0.01090

   et[1] here is a near-cancellation (ep ≈ −10.41, ev ≈ +10.43), which amplifies a ~5e-4 relative change in ep by ~1000×. 0.4.18 keeps a ~30% offset even at etol = 1e-12: the wrong super-diagonal is a systematic bias, not something that converges away.

   DIII-D-like — the 1.78 + 0.53i is the auto grid, not the fix. On fixed grids v0.4.20 is benign:

   grid	0.4.18	0.4.20
   ldp 512	0.80073	0.80056
   ldp 1024	0.80093	0.80062

   On the auto grid, 0.4.18 + only the fix and real 0.4.20 take the identical 822-step path, yet give 0.8006 vs 1.78 + 0.53i — a roundoff-level difference (e.g. the rewritten cubic eval kernel) producing an O(1) change with a large imaginary part. That's the same amplification family as #347, now in the forward-integrator energy itself, and it's what #349/#370/#525 describe: the fix changes the curvature measure, the auto grid gets coarser (598 → 320 surfaces, as noted above), and the forward solution on that grid is ill-conditioned. Downstream DIII-D quantities (resonant field, NTV torque ×1700, island widths) inherit it.

   Proposal: take the bump (it's a correctness fix), but sequence it with #525's knot-density floor so diiid_n1 doesn't re-pin the harness to a pathological value. Separately, a large Im(et[1])/Re(et[1]) seems worth a warning as a cheap detector for this failure mode.

   👍React with 👍1Reacted by Min-Gu Yoo
9. mgyoo86 commented on Oct 16, 2026

   

   mgyoo86

   on Oct 16, 2026

   Collaborator

   More actions

   Thanks @matt-pharr for the much more thorough investigation, this is very helpful.

   I mostly looked at the base quantities, which is why I estimated at most a few percent, and I didn't realize how much the near-cancelled values matter here. Either way, the fix is clearly needed.

   More generally, a few notes on sensitive systems like this, where ULP-level differences can get amplified:

   • Fundamentally, ULP-level results can't be trusted. Even the order of a simple sum changes the answer:

     0.3 - 0.4 + 0.1   # -2.7755575615628914e-17
     0.1 + 0.3 - 0.4   #  0.0
     0.1 - 0.4 + 0.3   # -5.551115123125783e-17

     So ULP-level results can change at any time, even in a minor or patch release of a package (FastInterpolations included). Simply reordering a single operation is enough.
   • Even the exact same source code can give different ULP-level results on different hardware, OS or Julia versions, depending for example on whether the compiler fuses a*b - c into an FMA (0.1*3.0 - 0.3 gives 5.55e-17, fma(0.1, 3.0, -0.3) gives 2.78e-17) or how SIMD reorders a reduction. For example, running the same DIII-D-like example with FastInterpolations of v0.4.20 on my machine (macOS arm64, Julia 1.13.2) gives 0.80062 rather than 1.78 + 0.53i. In practice, identical behavior at the ULP level is not something we can expect.
   • So any ULP-level analysis needs care, because results can differ from one person or environment to another. (For what it's worth, at that level v0.4.20 is slightly more exact: see fix(cubic): exact node values from the cubic value kernel ProjectTorreyPines/FastInterpolations.jl#211.)
   • Fundamentally, I think the direction should be to remove this kind of error amplification itself. A high-precision reference (e.g. BigFloat) would also help separate real errors from round-off amplification. FastInterpolations already supports BigFloat, but GPEC is largely hard-coded to Float64, so it would need a substantial refactor. Happy to help if there's interest.