Repository navigation
numr 0.10.0
A public numr::distance API on plain slices. SIMD kernels back it, and CPU cdist and pdist now use them.
Added
numr::distance— a public API on plain slices. It needs noTensororRuntime.- Float functions:
dot,l2_squared,manhattan, andcosine_distance. Each has an_f32and an_f64twin, such asdot_f32andcosine_distance_f64. - One query against many rows: the
*_many_*twins read a contiguous row-major buffer, such asdot_many_f32(query, rows, d, out). - Int8:
dot_i8,dot_i8_scaled, anddot_i8_many. The sum is exact and a total outside the i32 range saturates once at the end. - A length mismatch panics. The message names the function and both lengths.
- Cosine distance of a zero vector returns 0. For unit-length vectors, cosine distance equals
1 - dot.
- Float functions:
Kernelshandle —Kernels::detect(),Kernels::scalar(), andKernels::with_level(level) -> Option<Kernels>.with_levelreturnsNonewhen the CPU does not support the level. A caller can hoist dispatch or force a level withoutunsafe.numr::distance::SimdLevelis re-exported. It is#[non_exhaustive]and has no ordering, so a later release can add a level without a breaking change.- SIMD distance kernels — f32 and f64
dot, squared Euclidean, Manhattan, and cosine. Each has an AVX-512 path with a masked tail, an AVX2+FMA path, a NEON path, and a scalar fallback. NEON also serves NEON+FP16. Dispatch runs once per call on the cachedSimdLevel. Each ISA keeps several independent accumulators to hide FMA latency. benches/distance.rs— acdistbenchmark for squared Euclidean, cosine, and Manhattan at d = 128, 384, 768, and 1536.
Changed
- Minimum Rust version — 1.95, up from 1.89. A build with
--features f16on aarch64 does not compile on 1.89. The aarch64f16conversion kernels use NEON half-precision intrinsics that are stable only from 1.94. - CPU
cdistandpdist— Euclidean, squared Euclidean, Manhattan, and cosine route through the new kernels for F32 and F64.- They run in parallel across fixed work units. The unit size depends only on
d, so the output bits are identical for every thread count. - A single query against many rows (
n = 1) also runs in parallel. - F16, BF16, and FP8 inputs convert to F32 once, run on the F32 path, and convert back once.
- They run in parallel across fixed work units. The unit size depends only on
- Numerics — the summation order of these metrics differs from the old sequential loop and differs between SIMD levels. Results can differ in the last bits.
- Cosine keeps the
denom == 0 -> 0rule. - The CPU accumulator width still matches CUDA: f32 for narrow floats and f64 for F64. The CPU no longer sums term for term in the same order as CUDA.
- The CUDA and WebGPU parity tests for distance pass against the CPU path, on CUDA and on Vulkan through WGPU.
- Cosine keeps the
- Internal — block and row kernels replace the serial crate-internal
cdist_kernelandpdist_kernel. The public API does not change.
Removed
- cuBLAS —
CudaClientno longer creates a cuBLAS handle, so thecudafeature no longer loadslibcublasat run time. No numr kernel used it.CudaClient::cublas()andCudaError::CublasErrorare removed.- cudarc's default features are off. numr enables
driverandnvrtconly, so cuBLASLt, cuRAND and the CUDA runtime library are not built either.
Fixed
- Cosine and correlation distance at extreme magnitudes — the denominator was
sqrt(a * b)of two squared norms, which overflowed f32 once|a| * |b|passed about 1.8e19 and flushed to zero below about 3.7e-23.cosine_distance_f32(&[1e10], &[1e10])returned 1 instead of 0, and[1e-12]against[-1e-12]returned 0 instead of 2. The denominator is nowsqrt(a) * sqrt(b). The CPU, CUDA, and WebGPU kernels all changed, forcdist,pdist, andnumr::distance. Each squared norm must still fit the float type, so a component above about 1e19 (f32) or 1e154 (f64) still gives NaN. Results in the normal range differ from before only in the last bits. - Shape overflow — a tensor shape whose element count overflows
usizereturned a wrapped, too-small allocation.pdiston a[2^32 + 1, 0]tensor then wrote past the end of the output.Tensor::empty,Tensor::from_slice, and the scalar and typed fill constructors now returnError::InvalidArgumentwhen the count overflows.pdistandsquareform_inversedo the same for the pair count. - x86-64 SIMD detection —
SimdLevel::Avx512now also requires AVX2. A virtual machine that hides AVX2 but reports AVX-512 no longer reaches AVX2 code. - Build without
rayon—--no-default-featuresno longer warns about two matmul items that only the rayon path uses. - aarch64 clippy —
erf_f64usesFRAC_2_SQRT_PIinstead of a literal.
Performance
User-space instructions per cdist call, measured with perf stat -e instructions:u. Setup: n = 1, m = 64, one thread, AVX2+FMA host. Each count includes about 3,000 instructions of tensor allocation and dispatch. These are instruction counts, not wall-clock times.
| Metric | d | Before | After |
|---|---|---|---|
| Squared Euclidean | 128 | 43,453 | 10,625 |
| Squared Euclidean | 384 | 121,265 | 19,317 |
| Squared Euclidean | 768 | 238,011 | 34,506 |
| Squared Euclidean | 1536 | 471,485 | 61,168 |
| Cosine | 128 | 103,473 | 15,464 |
| Cosine | 384 | 300,093 | 29,033 |
| Cosine | 768 | 595,005 | 51,721 |
| Cosine | 1536 | 1,184,828 | 93,534 |
| Manhattan | 128 | 53,628 | 11,561 |
| Manhattan | 384 | 151,932 | 22,313 |
| Manhattan | 768 | 299,388 | 41,088 |
| Manhattan | 1536 | 594,301 | 73,121 |
Known limits
- The AVX-512 kernels compile and pass clippy. They have never run. The development host has no AVX-512, and QEMU TCG does not emulate it.
- A separate read of the code traced the AVX-512 loop bounds, tail masks, and feature requirements by hand.
- The NEON kernels ran only under
qemu-aarch64(user-mode emulation). They have not run on real ARM hardware. - Instruction counts for NEON and AVX-512 are not measured.
Full diff: v0.9.0...v0.10.0