compiler: compute sparse position/floor in fp64 to fix off-by-one cell shift - #2992
compiler: compute sparse position/floor in fp64 to fix off-by-one cell shift#2992mloubout wants to merge 1 commit into
Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #2992 +/- ##
==========================================
- Coverage 83.54% 78.23% -5.31%
==========================================
Files 257 251 -6
Lines 53922 53432 -490
Branches 4613 4600 -13
==========================================
- Hits 45047 41803 -3244
- Misses 8076 10764 +2688
- Partials 799 865 +66
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:
|
24006cd to
04411cd
Compare
|
Check out this pull request on See visual diffs & provide feedback on Jupyter Notebooks. Powered by ReviewNB |
1ad6fd3 to
038d26b
Compare
FabioLuporini
left a comment
FabioLuporini
left a comment
There was a problem hiding this comment.
I think the main question that I have is: why not inlining DOUBLE(...) casts directly when building the interpolation operation? why is this implemented via hoisting into fp64 scalars?
| yield self._do_generate(exprs, exclude, cbk_search) | ||
|
|
||
|
|
||
| def _is_floor(e): |
There was a problem hiding this comment.
this is quite hacky, let the search look for the specific type if it really needs to be that specific
| subs[inner] = DOUBLE(v) | ||
|
|
||
| exprs = [uxreplace(e, subs) for e in exprs] | ||
| return exprs, aliases |
|
|
||
| # Predicate on Cluster used to pick which ones this pass fires on. | ||
| # Subclasses override to target a different kind of cluster. | ||
| _cluster_filter = staticmethod(lambda c: c.is_dense) |
There was a problem hiding this comment.
why don't we always do it with all floors (and maybe not just that), be it dense or sparse ?
| issubclass(k.dtype, np.integer)): | ||
| continue | ||
| inner = k.base | ||
| if _is_floor(inner) and (v.free_symbols & aliaseds): |
There was a problem hiding this comment.
instead of just floors, we could name this class CireInvariantGeometry or something like that, and have it target all step-wise math operations that are so accuracy-sensitive
you could have for example a _types class attribute, and use it here
038d26b to
307d39c
Compare
The linear interpolator computes `(c - o)/h` in the grid's dtype, which for fp32 grids rounds coord/origin/spacing to fp32. That rounding can push the position across an integer boundary, so `floor((c - o)/h)` picks a different cell than the fp64-truth cell -- and CPU vs GPU can disagree on the fractional part while agreeing on the integer position, producing inconsistent injection weights. Compute `pos = (c - o)/h` and `floor(pos)` in fp64 by casting the fp32 free symbols and by substituting the spacing symbols with their fp64 decimal value (recovered via the fp32 short-form round-trip). A new CIRE subpass `CireInvariantsSparse` hoists `pos` and `floor(pos)` out of the per-stencil-point inner loop into per-source preamble Arrays. The floor tab is stored as int32 (half the memory of an fp64 tab, no precision loss), and bare `floor(pos)` uses inside `pos - floor(pos)` are rewritten to `DOUBLE(int_tab)` so both consumers share the tab. The `sympy_dtype` inference is taught to recognize `Cast` (outermost and inner) so the printer emits `floor` (fp64) instead of `floorf` (fp32) when a `DOUBLE(...)` cast is present in the expression.
307d39c to
f25a1eb
Compare
2 participants
The lifiting is implemented as a Cire pass this way it'll catch generic floor/... as well