FD weights: 17-digit C literals + canonicalized Taylor weights (fixes #2976)#2977
Open
ggorman wants to merge 1 commit into
Open
FD weights: 17-digit C literals + canonicalized Taylor weights (fixes #2976)#2977ggorman wants to merge 1 commit into
ggorman wants to merge 1 commit into
Conversation
…2976) _PRECISION 9 -> 17: 17 significant digits is the minimum that round-trips an IEEE-754 double, so the C literals parse back to the exact rational weight (sub-ULP) instead of carrying a ~1.7e-10 relative truncation that silently floors fp64 cross-implementation checks. Raising the precision exposed a latent route-inconsistency the 9-digit rendering had been masking: Taylor weights built through a float x0 shift (e.g. div(f, shift=.5)) accumulate ~1e-14 relative float error, so the same stencil built through two routes compared unequal at 17 digits. Taylor-computed weights are now canonicalized (nearby-rational recovery within 5e-14 relative, then evalf) so identical stencils are identical expressions; user-supplied weights are rendered as-is. Tests: test_derivatives now imports the library _PRECISION instead of a local copy; 9-digit string expectations regenerated at 17 digits; new TestCoefficientPrecision regression class (weight-literal fp64 round-trip + evaluated-derivative Float atoms round-trip to their exact rationals).
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #2977 +/- ##
==========================================
- Coverage 83.54% 83.48% -0.06%
==========================================
Files 257 257
Lines 53591 53623 +32
Branches 4585 4588 +3
==========================================
- Hits 44770 44766 -4
- Misses 8028 8053 +25
- Partials 793 804 +11
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
mloubout
reviewed
Jul 21, 2026
| r = sympy.Rational(a).limit_denominator(10**12) | ||
| except (TypeError, ValueError, ZeroDivisionError): | ||
| return a | ||
| if a == 0 or abs(float(r) - float(a)) <= 5e-14 * abs(float(a)): |
Contributor
There was a problem hiding this comment.
float isn't safe its value will depend on the compilation flags and platform of python. This should likely use np.isclose with proper atol/ftol and np.float64 instead of float
| at ``_PRECISION`` so the emitted C literal round-trips the intended double | ||
| exactly. | ||
| """ | ||
| def _canon(a): |
Contributor
There was a problem hiding this comment.
why is this a local function?
| non-rational user-supplied coefficients are value-preserved), then render | ||
| at ``_PRECISION`` so the emitted C literal round-trips the intended double | ||
| exactly. | ||
| """ |
Contributor
There was a problem hiding this comment.
this is quite unreadable honestly, no idea what it's trying to say.
| # compare unequal; user-supplied weights are rendered as-is. | ||
| scale = dim.spacing**(-deriv_order) if scale else 1 | ||
| weights = [sympify(scale * w).evalf(_PRECISION) for w in weights] | ||
| if computed_weights: |
Contributor
There was a problem hiding this comment.
This flag is quite ugly, would make more sense to just have the _canonicalize_weight directly part of the weight computation. This "coefficients == 'taylor'" also ignore custom rules and things like the drp mode
| from devito.warnings import DevitoWarning | ||
|
|
||
| _PRECISION = 9 | ||
| from devito.finite_differences.finite_difference import _PRECISION # noqa |
Contributor
There was a problem hiding this comment.
No blank line, wrong place (goes before devito.types)
| ]) | ||
| def test_fd_weight_literal_roundtrips_fp64(self, p, q): | ||
| w = sympy.Rational(p, q) | ||
| emitted = float(sympy.sympify(w).evalf(_PRECISION)) |
| f"{p}/{q}: evalf(_PRECISION={_PRECISION}) -> {emitted!r} != {exact!r} " | ||
| f"(rel err {abs(emitted - exact) / abs(exact):.3e})") | ||
|
|
||
| def test_evaluated_derivative_floats_roundtrip_fp64(self): |
Contributor
There was a problem hiding this comment.
Not sure what this tests, or why it's testing a roundtrip for fp64 with an fp32 grid/....
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Fixes #2976.
What
_PRECISION9 → 17 (finite_differences/finite_difference.py). 17 significant digits is the minimum that round-trips an IEEE-754 double, so the C literal parses back to the exact rational weight (sub-ULP) instead of carrying a ~1.7e-10 relative truncation. No runtime cost — only the codegen text changes.Taylor-weight canonicalization. Raising the precision exposed a latent inconsistency the 9-digit rendering had been masking (this is what the old "avoid roundup errors" comment was living with): Taylor weights built through a float
x0shift (e.g.div(f, shift=.5)) accumulate ~1e-14 relative float error in the weight recursion, so the same stencil built through two routes compared unequal at 17 digits (0.041666666666666664vs…667, caught bytest_shifted_div). Taylor-computed weights are now canonicalized — nearby-rational recovery within 5e-14 relative viaRational(w).limit_denominator(1e12), thenevalf(_PRECISION)— so identical stencils are identical expressions and every emitted literal is the correctly-rounded double of the exact rational. User-supplied weights (weights=/ symbolic coefficients) are rendered exactly as before — canonicalization applies only to the computed Taylor path.Tests.
test_derivatives.pynow imports the library_PRECISIONrather than pinning a local copy (so these exact-equality tests track the library value, liketest_symbolic_coefficients.pyalready does); the four 9-digit string expectations are regenerated at 17 digits; and a newTestCoefficientPrecisionclass pins the invariant: every FD weight literal round-trips to its exact rational as a double (this class fails at_PRECISION = 9— e.g.8/315→2.53968254e-2, 1.7e-10 off — and passes at 17).Why it matters
In a downstream verification corpus we swept 332 Devito-vs-reference gates and found 11 (in 9 projects) whose tolerances had been silently loosened into the 1e-7..1e-9 range specifically to swallow this floor, several with in-code comments like "Devito's generated C rounds the rational stencil coefficients". Each such gate is blind to real coefficient defects at the magnitude its tolerance implies. Details and measurements in #2976.
Test status
Ran locally (macOS, fp64):
tests/test_derivatives.py+tests/test_symbolic_coefficients.py— 409 passed;tests/test_unexpansion.py25 passed;tests/test_caching.py60 passed. A repo-wide grep found no other 9-digit literal expectations. Please run the full CI matrix — happy to iterate on anything it surfaces.Note: to be reviewed and merged by the devito maintainers (@mloubout) — not by us.