Repository navigation
Native numpy: float32 and int32 arrays - #228
Merged
Merged
Conversation
Python files that use numpy's core now build with a numpy written in Rust (src/dynpy/rt_np.rs, rt_npfn.rs, rt_npfmt.rs): float64/int64/bool arrays with NEP 50 promotion, strided views, masks and integer-array indexing, broadcasting, operators, ~60 ufuncs, reductions in numpy's pairwise order, sort/unique/stacking, and numpy 2.x printing byte for byte. Anything it lacks (np.random, np.linalg, other dtypes) is refused at build time and the program gets CPython's numpy through the typed build, as before. The dynamic build no longer calls CPython's numpy itself (5 of 40 corpus programs gave wrong output that way). BREVI_NUMPY=0 turns it off. tests/numpy: 11 programs whose .expected files come from real numpy; the check fails if any of them runs off Brevi's numpy. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
src/dynpy/rt_nprand.rs: default_rng(seed) as numpy does it (SeedSequence, PCG64, the ziggurat normal and exponential with numpy's tables, Lemire bounded integers) and the legacy np.random.seed API (MT19937 init_genrand, polar gauss with its cached value, masked randint). random, integers, normal, uniform, standard_normal, exponential, choice (replace=False, p=), shuffle, permutation, rand, randn, randint, RandomState. Other distributions and ufunc methods (np.maximum.accumulate) are refused at build time, so the program gets CPython's numpy. Unseeded generators take OS entropy without file access. tests/numpy: 6 random programs, expected output from real numpy. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
…apezoid, ufunc reduce/accumulate Float vector np.dot now follows OpenBLAS's own summation order for the CPU (AVX2 or AVX-512, with its thread split), so it matches numpy to the last digit. det/solve/inv are written but stay off (go to CPython numpy) until they match OpenBLAS's LU bit for bit. Fixes str * np.int64 raising TypeError. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Follows OpenBLAS's single-thread LU step by step (getrf_single recursive blocking, getf2, dgemv_n, the SkylakeX/Haswell dgemm kernels' k-loop splits, trsm_kernel_LT/LN, trsv, daxpy), so results equal numpy's to the last digit on AVX-512 and AVX2 CPUs. slogdet returns SlogdetResult; matrix_power takes negative powers. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
poisson, binomial, geometric, negative_binomial, hypergeometric, zipf, logseries, gamma, standard_gamma, beta, chisquare, noncentral_chisquare, f, noncentral_f, standard_t, standard_cauchy, lognormal, laplace, logistic, gumbel, weibull, pareto, power, rayleigh, wald, vonmises, triangular, multinomial and dirichlet on Generator and legacy RandomState, ported from numpy's distributions.c / legacy-distributions.c with _common.pyx's argument checks. normal/uniform/exponential move onto the same parameter code (NaN scale now gives nan as in numpy). Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Models OpenBLAS dsyrk (single-thread and threaded drivers, LOWER kernel with 32x32 diagonal blocks) so A @ A.T, A.T @ A and np.dot forms match numpy to the last digit on AVX2 and AVX-512. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
Ports the LAPACK routines numpy's dgelsd runs (QR/LQ, dgebd2, dbdsqr, dlalsd's small path, dlartg/dlasv2/dlas2, dlascl) with the BLAS calls following OpenBLAS's kernels, including dnrm2's x87 sums (run with the x87 instructions on x86-64). polyfit supports full, w, rcond and cov; lstsq returns its rank as numpy.int32. Hot numpy kernels now use the CPU's FMA instruction instead of the software routine. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
numpy's own lstsq/polyfit digits differ between OpenBLAS's AVX-512 and AVX2 kernels, and the CI runner has no AVX-512. Brevi matches numpy on both; a CPU without AVX-512 now checks NAME.avx2.expected when one exists. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
float32/int32 dtypes in constructors, astype, arange, linspace and the np.float32/np.int32 scalars, with NEP 50 promotion, float32 arithmetic rounded once per operation, int32 wrapping with numpy's warnings, float32 reductions and float32 printing. @/dot, linalg, ufunc reduce/accumulate and SVML-backed float32 math stay with CPython's numpy (build gate plus a runtime refusal). Fixes found by fuzzing against numpy: Dragon4 ties in printed digits (and Python float repr ties), float32 exp-format cutoff, np.rint of ints, -0.0 sums, array ** 0.5/-1 as sqrt/reciprocal, shape= in big-array repr, int() of nan numpy floats, arange with a float step and int dtype. numpy's unstable default sort order of equal values is refused unless kind='stable'. CPython bridge: np.where and x ** -1. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_0136hkZgpjat3dn8BGZPJVF2
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.
Requested by Adam · project thread
What
PR 7a of native numpy:
float32andint32arrays and scalars.dtype=np.float32/np.int32(also"f4","i4",np.single,np.intc) work in constructors,astype,arangeandlinspace.np.float32(x)andnp.int32(x)work too.OverflowError.//and%use float32 divmod.np.float32(...)scalar repr, anddtype=in array repr.@/dot(sgemm),np.linalg, ufuncreduce/accumulate, and float32exp/log/trig/powers, where numpy uses SVML. A build gate sends these files to CPython's numpy, and a runtime refusal catches what reaches them anyway.Fixes found by fuzzing against numpy
240 random programs mixing float32, int32, float64 and int64 now match numpy byte for byte; the only differences are explicit refusals. Several of these bugs were already present for float64:
reprof floats in typed and dynamic builds (-254.62014770507812).np.rintof an int array now returns float64.-0.0values gives0.0(numpy's add.reduce starts from the identity).** 0.5/** -1/** 2with a scalar exponent runs numpy'ssqrt/reciprocal/square. This is numpy's fast_scalar_power, so(-0.0) ** 0.5is-0.0, and the warnings name the right ufunc.shape=(...).int()of a numpy float that is nan or inf now raises; a large value stays exact.np.arangewith a float step and an int dtype now fills as numpy does.MIN // -1, and-/absof the minimum int scalar, now warn.NPY_DISABLE_CPU_FEATURES). Withoutkind="stable",np.argsortof repeated values and sorting0.0with-0.0now raise NotImplementedError instead of silently giving a stable order.kind=is now accepted bynp.sort.np.wherewas lexed aswhere_, andx ** -1calledfloat()on the array.Tests
tests/numpy/narrow.py. Its output is the same with numpy's AVX2 dispatch and with OPENBLAS_CORETYPE=Haswell.tests/check.shpasses.sh tools/ci_local.sh: all checks passed (2003 baseline solutions, 0 regressions).Generated by Claude Code