Repository navigation
Preserve constant offsets in spline and cspline endpoint derivatives - #308
Conversation
The 'extrap' end derivatives were evaluated as SUM(c*f) with c(1) = -sum(c(2:4)). For a constant profile, or any profile with a large offset, on a grid clustered at an endpoint the four terms cancel catastrophically (error ~ eps*|f|/h^3). Evaluate the same four-point derivative from differences to the endpoint value.
logan-nc
left a comment
There was a problem hiding this comment.
Sounds good to me and I know @matt-pharr has issues with the ldp grid packing in the past. I haven't independently checked this though.
|
@matt-pharr I wanted to ping you on this since I know you've had ldp spline packing issues in the past. You are assignee just to make sure you see it. You can hit merge or do a more thorough review - challenge by choice. I was convinced by the included test but didn't independently check. |
|
@logan-nc it makes me very uncomfortable that we are changing the splines... I will need to look into this very carefully before merging |
|
@krystophny FYI a bunch of your recent PRs have been into master, they all should be into develop, and you should be using develop in general as that is the most up to date version. Master is for releases. |
Understood. I will make Claude understand that too. |
|
@logan-nc I had claude run some tests on the errors with the spline packing, and it seems like this hearkens way back to the classic spline investigations by @zrwang2015... #68 #55. Here are some test representations of different functions with right-hand-side packing:
It also found the reason that zhirui's splines were producing different/worse results in some cases, that being that it uses real types where it should use complex types in cspline_fit_classic and therefore throws away complex derivatives. So I am going to work that fix into this PR. I also checked out chris' changes and they do seem correct. I suggest merging this after I clean up the extraneous python tests. Funny enough this change would have saved me from the year of despair at the RDCON 'infinite surface shielding' issue :') |
… be appropriately complex
…endpoint-cancellation
matt-pharr
left a comment
There was a problem hiding this comment.
Removed the tests and added the type fix for the classic splines. This is merge-able now.
|
Thanks @krystophny for this deep fix and @matt-pharr for the detailed checking / tracking. This is definitely great to have, and hopefully no more researchers lose a year to faulty spline! The classic spline fix is great to have too. But those ate definitely slower, so this one is a value add beyond that 😀 |
Fixes #306.
spline_fit_ahgandcspline_fit_ahg(endmode "extrap") computed the end derivative asSUM(c*f)withc(1) = -sum(c(2:4)). The terms are O(f/h) and cancel, so a constant or large-offset profile on a grid clustered atthe right end acquires an error of order eps*|f|/h (4.8e-2 for h ~ 1e-5, |f| ~ 1). The change evaluates the same
four-point derivative from differences to the end value, in
equil/spline.fandequil/cspline.f. No otherbehaviour changes.
Test:
regression/test/test_spline_endpoint_polynomials.pyfits constant and cubic profiles (and a complex constantprofile with offset 300) on the 'ldp' sine-squared grid with 512 and 1024 intervals and checks the end derivatives,
the end values and the integral through
spline_int. The test compiles onlyio,local,utils,splineandcsplinefromequil/. Before the change it fails at 512 and 1024 intervals for every profile (largest derivativeerrors: 0.11 and 12.4 for the cubic, 0.89 for the constant, 586 for the complex constant); after the change it passes.
End to end (lar deck of the issue, A = 40, 3/1 drive): DCON plasma energy and GPEC response are the same for
mpsi = 256, 320, 384, 512 (energy -6.434e-2, difference below 1e-3 relative) and for ma = 512, 1024; for
A = 70..1000 the total free-boundary energy is positive and scales as 1/A.
Chris&AI