Repository navigation
Write expr[:] into the NumPy result for N-D arrays; reuse chunk buffers in chunked paths - #753
Merged
FrancescAlted merged 4 commits intoNov 3, 2026
Conversation
Two common cases fell off miniexpr, and without numexpr they now run chunk by chunk through evaluate_chunk(): - Full std()/var() of float32/float64 (no axis, no where, same dtype) reduced (x - mean) ** 2 with the mean as a broadcast (1, ..., 1) operand. They are now two single-pass reductions: the mean, then the mean of (x - m) * (x - m) with m inlined as a scalar. The sign of m is folded into the operator, as miniexpr reductions mis-evaluate a negative literal on the right of a binary operator. Other cases keep the general path; std() and var() share one implementation. - where(cond, x, y) with a scalar x or y: specialize_miniexpr_inputs() inlined the scalar and dropped its operand before the where() was fused in, so the fused expression referenced a missing operand, miniexpr failed to compile it and the evaluation fell back. The where() is now fused before specializing. - The bare-array LazyExpr that the NDArray reductions build takes the array's dtype instead of inferring it. Measured on a Ryzen 9 7950X (300x512x512 float32 in memory, 28 threads), against this branch at 7db7b09: a.std() 470 -> 53 ms, (a.mean(), a.std()) 498 -> 78 ms, where(a > 100, a, 0).compute() 690 -> 469 ms; 100x100 a.std() 0.73 -> 0.52 ms. std/var may differ from before in the last ulps. Taken from the lazy-expression commit of Blosc#744 (b28fb0e), which was deferred until this branch lands. Co-Authored-By: Claude Opus 5.5 <[email protected]>
Member
|
This LGTM, but will need to wait until #752 would be merged. |
Fusing where(cond, x, y) with scalar x/y (previous commit) exposed miniexpr to cases it gets wrong, as Francesc found in review: - where(a > 0.1, 1, 0) with a float32 `a` of 0.1 selected 1. miniexpr types the literals of the fused expression by its output dtype and computes in it: with an int64 output 0.1 became 0, with a float32 output it compared against the double 0.1, where NumPy rounds it to float32 first. Arithmetic in the condition also ran in the output dtype (where(a + b == c, 1, 0) with float32 operands). - where(a > 5, np.uint64(2**63 + 1), a) gave 2**63: miniexpr parses literals into C doubles. Inside where() it also passes 64-bit integer operands and branch values through doubles. where() is now fused only when computing the condition in the output dtype changes nothing: all its operands have the output dtype, or it only compares operands that cast safely to it against literals their dtypes hold exactly; never with 64-bit integers. Otherwise it falls back to the chunked evaluation. This also fixes array branches, which were fused before and had the same problems (where(a > 0.1, a, b)). Masks, label maps, thresholds and same-dtype arithmetic stay on miniexpr; where(a > 100, a, 0).compute() keeps its speed. Expressions with integer literals beyond 2**53, written or inlined scalars, are no longer sent to miniexpr (also a + np.uint64(2**63 + 1)). Tests: the review examples, a NumPy differential grid over condition dtypes, literals, comparisons, output dtypes and scalar/array branches, arithmetic conditions, large integer constants, and checks that the common mixed-dtype where() patterns still compile with miniexpr. Co-Authored-By: Claude Opus 5.5 <[email protected]>
FabianIsensee
force-pushed
the
fast-lazy-getitem-reductions
branch
from
November 3, 2026 13:25
5ca92ec to
dbff1c9
Compare
Two CI failures of the previous commits: - NumPy 1.26: a.std(ddof=1) and a.var(ddof=1) of float32 returned float64. The full std/var returns a float32 scalar, and under NumPy 1.x's value-based casting a float32 scalar times the Python int count of the n / (n - ddof) rescaling becomes float64 (with keepdims, a float32 array times a count above 65535 does too). The rescaled variance is now cast back to its own dtype; with NumPy 2 this changes nothing. This fixes the 24 test_std_var_full_reduction[...-1-...- float32] failures. - WebAssembly: test_large_integer_constants_exact and test_where_literal_condition_matches_numpy added 2**53 to a plain np.arange(), which is int32 on wasm32 and raised OverflowError. They build int64 inputs explicitly now. Not changed here: on NumPy 1.26, mean() of float32 already returns float64 for full reductions and for axis reductions over more than 65535 elements, and mean(where=) does on NumPy 2 too (also on 7db7b09). Co-Authored-By: Claude Opus 5.5 <[email protected]>
…rs in chunked paths - fast_eval()'s NumPy-output path wrote the miniexpr results of each chunk straight into the result only for 1-D arrays. It now does so for any array whose chunks are consecutive runs of the C-order result: blocks C-contiguous in their chunk, chunks C-contiguous in the array, only the last one partial (_chunks_tile_flat()), as with the default partitions. An N-D `out` must be C-contiguous for that; a 1-D one may still be strided. Other arrays keep the uncompressed temporary NDArray. - Chunked evaluation and reductions (axis=..., argmin, cumulative_sum, ...) decompress the chunks of operands whose blocks are C-contiguous runs of their chunk straight into reused chunk-sized buffers (_chunk_buffers()), instead of allocating per chunk; the unused `aligned` plumbing is gone. blocks_contiguous() moves to module level in ndarray.py for that. - One compiled regex for all constructor checks. Measured on a Ryzen 9 7950X (300x512x512 float32 in memory, 28 threads), on top of the std/var and where() change: (a + b)[:] 116 -> 71 ms, a.sum(axis=0) 29.3 -> 26.3 ms. Taken from the lazy-expression commit of Blosc#744 (b28fb0e), which was deferred until this branch lands; its own NumPy-output mechanism is replaced by extending the one of this branch. Co-Authored-By: Claude Opus 5.5 <[email protected]>
FabianIsensee
force-pushed
the
fast-lazy-getitem-reductions
branch
from
November 3, 2026 15:12
dbff1c9 to
b08bfd9
Compare
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
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.
Hi Francesc,
this is the second PR into #741: the rest of the lazy-expression speedups of #744.
Depends on #752: it is stacked on it, so until that is merged, this PR also shows its commit. The new commit is the last one.
What:
expr[:]:fast_eval()'s NumPy-output path wrote each chunk's miniexpr results straight into the result only for 1-D arrays. It now does so for any array whose chunks are consecutive runs of the C-order result: blocks C-contiguous in their chunk, chunks C-contiguous in the array, only the last chunk partial (_chunks_tile_flat()), which is what the default partitions give.outmust be C-contiguous for that; a 1-D one may still be strided, as before.axis=...,argmin,cumulative_sum, ...): operands whose blocks are C-contiguous runs of their chunk now decompress straight into reused chunk-sized buffers, instead of allocating per chunk. The unusedalignedplumbing is gone, andblocks_contiguous()moves to module level inndarray.py.Measured on a Ryzen 9 7950X: 300x512x512 float32 in memory, 28 threads, against c-blosc2
main; medians of two interleaved rounds.(a + b)[:]a.sum(axis=0)a.std()where(a > 100, a, 0).compute()Tests:
test_getitem_into_numpy_layouts: layouts that do and do not tile, partial chunks, float32/int16/complex128, and that the result lands in the direct buffer exactly when the layout tiles.test_reduction_axis_chunk_layouts:sum/max/argmin/cumulative_sumover several axes and chunk/block layouts.test_evaluate_numpy_prefilter.pypasses, including the strided 1-Doutcases.pytest tests: 10788 passed. The 2 failures are theb2viewhttp tests that need aiohttp, as onmain.AI disclosure: the port, the tests and the measurements were done with Claude Code; I reviewed them.
Best,
Fabian
🤖 Generated with Claude Code