diff --git a/CHANGELOG.md b/CHANGELOG.md index b2a9209..d6000c7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -13,6 +13,66 @@ against, and semantic-release manages only the **patch** component. ## [Unreleased] +### Added + +- **MTL5 pinned to v5.11.0, and `accumulator='i32'`** + ([#88](https://github.com/stillwater-sc/mtl5-python/issues/88) phase 2). The + narrow integer element types phase 1 registered can now compute: + + ```python + a = np.random.randint(0, 256, 4096, dtype=np.uint8) # activations + w = np.random.randint(-128, 128, 4096, dtype=np.int8) # weights + mtl5.mixed.dot(a, w, accumulator="i32") + ``` + + MTL5 v5.11.0 routes `dot` over 8- and 16-bit operands onto the + hardware widening multiply-accumulate — `vpmaddwd` / `vpdpbusd` on x86, + `SMLAL` / `SDOT` on NEON — and an int32 accumulator is what those + instructions accumulate into. + + **Every signedness pairing is accepted, in either order.** `u8 × i8` is VNNI's + native shape on x86 and what quantized inference is written in; ARM implements + the symmetric pairings first. A dot product is symmetric, so MTL5 swaps the + operands onto whichever form the machine has. The kernel below + (`simd::reduce_dot_widen`) rejects `(int8, uint8)`, but that restriction is the + kernel's and re-exposing it would refuse a call the library can serve. + + **The overflow contract is part of the API, not a footnote.** Products are + always exact; the sum wraps, and how soon depends on operand *magnitude* + rather than vector length — roughly `2^(31-2b)` terms at `b` bits. Measured at + full range: one `i16 × i16` product uses 2³⁰ of the int32 range, so **two** + already overflow it, while `i8 × i8` holds **131071** terms. That five order + of magnitude gap is why the quantized-inference instructions are 8-bit. The + wrap is two's complement and therefore bit-identical across lane counts, + backends and thread partitions — reproducible, but still wrapping. It is + stated in the docstring and asserted in the tests rather than left to be + discovered. + + Two guardrails, both deliberate. `accumulator=` is **required** for these + dtypes, unlike every other: omitting it means element precision, which is + exactly the wrapping phase 1 refused to expose, and silently redefining `None` + for three dtypes alone would be worse than asking. And `result='element'` is + **refused**, since rounding an int32 sum back to an 8-bit element + re-introduces the wrap the accumulator exists to avoid. + + `norm` and `frobenius_norm` remain unregistered for these types: `two_norm` + takes `sqrt` of the accumulated sum, and the API has no way to say + "accumulate in int32, deliver a real square root". + + `mtl5.dtypes()` is deliberately **unchanged**. Its contract is the set + `convert()` can target — the suite parametrizes over it and converts into + every entry — and `convert()` cannot target an integer: re-quantizing reals + into 8 bits is a quantization scheme (scale, zero point, rounding mode), not + a cast, and a naive version would silently clip everything outside + [-128, 127]. `mtl5.mixed.accumulators('i8')` answers for them instead. + + The v5.11.0 pin is what makes this sound rather than merely available. + `batch` does not exist before it, so the kernels are absent — and + `detail/wrapping_arithmetic.hpp` is new in it, without which the generic + integer loops are **UB** on overflow rather than the documented modular wrap. + For these operand widths overflow is the normal regime, not an edge case. + + ### Added - **Narrow integer element types: `i8`, `i16`, `u8`** diff --git a/CMakeLists.txt b/CMakeLists.txt index 489eb41..65a3606 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -100,17 +100,26 @@ endif() # regression on the refactor path, with no diagnostic. # v5.10.0 -- mtl5_module.cpp calls mtl::util::build_isa_list() (MTL5 #443), # which does not exist before this tag. +# v5.11.0 -- mtl5_mixed_precision.cpp offers accumulator='i32', which is the +# widening integer dot: int32 accumulation over 8- and 16-bit +# operands, on vpmaddwd / vpdpbusd / SDOT. Two things arrive with +# this tag and neither is optional. batch does not exist +# before it, so the kernels are simply absent. And +# detail/wrapping_arithmetic.hpp is new here: without it the +# generic integer loops are UB on overflow rather than the +# documented two's-complement wrap, and overflow is the NORMAL +# regime for these operand widths, not an edge case. # # The version argument to find_package is what applies that floor to the OTHER # build path. A bare find_package(MTL5 QUIET) accepts *any* system-installed # MTL5 and skips the FetchContent block entirely, so a developer with 5.7.x # installed would silently build against it and get ~140 template errors deep # inside norms.hpp instead of a version message. MTL5 ships its ConfigVersion -# with COMPATIBILITY SameMajorVersion, under which 5.10.0 rejects an installed -# 5.9.x and accepts 5.10.0 or 5.11.0 -- exactly the floor we want. (It also +# with COMPATIBILITY SameMajorVersion, under which 5.11.0 rejects an installed +# 5.10.x and accepts 5.11.0 or 5.12.0 -- exactly the floor we want. (It also # rejects a 6.x install, which is correct: an MTL5 major bump is the manual # intervention case in the version policy, not something to absorb silently.) -find_package(MTL5 5.10.0 QUIET) +find_package(MTL5 5.11.0 QUIET) if(NOT MTL5_FOUND) include(FetchContent) # Suppress MTL5's own tests, examples, and install targets @@ -120,7 +129,7 @@ if(NOT MTL5_FOUND) FetchContent_Declare( mtl5 GIT_REPOSITORY https://github.com/stillwater-sc/mtl5.git - GIT_TAG v5.10.0 + GIT_TAG v5.11.0 GIT_SHALLOW TRUE EXCLUDE_FROM_ALL ) diff --git a/README.md b/README.md index 430a3e2..cdbb9f4 100644 --- a/README.md +++ b/README.md @@ -115,6 +115,39 @@ quires have known upstream limitations documented in `accumulator=` is available on `dot`, `norm` (ord=2), `frobenius_norm`, `matvec` and `matmul`. +### Integer operands: `accumulator="i32"` + +The same idea with integers, which is where the hardware is. `int8` and `int16` +operands accumulated in `int32` map onto the widening multiply-accumulate — +`vpmaddwd` / `vpdpbusd` on x86, `SMLAL` / `SDOT` on NEON: + +```python +a = np.random.randint(0, 256, 4096, dtype=np.uint8) # activations +w = np.random.randint(-128, 128, 4096, dtype=np.int8) # weights +mtl5.mixed.dot(a, w, accumulator="i32") +``` + +`u8 × i8` is VNNI's native pairing on x86; ARM implements the symmetric ones +first. Either order works — a dot product is symmetric, so MTL5 swaps the +operands onto whatever the machine has. + +**The sum wraps, and sooner than vector length suggests.** Products are always +exact, but headroom goes as operand *magnitude*: about `2^(31-2b)` terms at `b` +bits. Measured at full range, one `int16 × int16` product uses 2³⁰ of the int32 +range so **two** overflow it, while `int8 × int8` holds **131071** terms. That +gap is why quantized inference is 8-bit. The wrap is two's complement, hence +bit-identical across lane counts, backends and thread counts. + +`accumulator="i32"` is required for these dtypes — the default is element +precision, and an 8-bit accumulator overflows almost immediately. Build the +arrays with NumPy: `convert()` does not target integers, because re-quantizing +reals into 8 bits is a quantization scheme rather than a cast. + +A released wheel is built at the x86-64 baseline, so it gets int8's **bandwidth** +win — one byte per element against float64's eight — but not the VNNI +instruction. `mtl5.build_info()["build_isa"]` reports which you have; build with +`-C cmake.define.MTL5_NATIVE_ARCH=ON` to reach it. + ### Iterative refinement Factor cheaply in a low precision, then recover accuracy with a residual formed diff --git a/pyproject.toml b/pyproject.toml index 8dff376..1f56a58 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -47,7 +47,7 @@ build-backend = "scikit_build_core.build" [project] name = "mtl5" -version = "5.10.3" +version = "5.11.0" description = "Python bindings for MTL5 — NumPy/SciPy/JAX/PyTorch interop with hardware accelerator dispatch" readme = "README.md" license = {text = "MIT"} diff --git a/python/src/mtl5_mixed_precision.cpp b/python/src/mtl5_mixed_precision.cpp index a654fd4..344cb07 100644 --- a/python/src/mtl5_mixed_precision.cpp +++ b/python/src/mtl5_mixed_precision.cpp @@ -127,6 +127,20 @@ double dispatch_acc(AccKind kind, F&& f) { return f.template operator()::type>(); else throw std::invalid_argument("accumulator='quire' unavailable for this dtype"); + case AccKind::I32: + // Guarded exactly as Quire is, and for the same reason: a discarded + // if-constexpr branch is not instantiated, so `mtl::dot` + // over a posit or cfloat vector is never formed. parse_acc rejects + // 'i32' for those dtypes first; this is the second line of defence + // and the one the compiler enforces. + // + // int32 is exactly representable in a double's 53-bit significand, + // so returning it through this function's double is lossless. That + // is worth stating because it looks like a narrowing bug and is not. + if constexpr (std::is_integral_v) + return f.template operator()(); + else + throw std::invalid_argument("accumulator='i32' unavailable for this dtype"); } throw std::invalid_argument("unreachable accumulator kind"); } @@ -314,6 +328,169 @@ void register_mixed_native(nb::module_& mx) { }, "A"_a, "accumulator"_a = nb::none()); } +/// Dot over two DIFFERENT narrow integer operand types. +/// +/// The signedness pairing is not a detail: `u8 x i8` -- unsigned activations +/// against signed weights -- is VNNI's NATIVE shape on x86 and the one quantized +/// inference is written in. A same-type-only surface would put the whole point +/// of the 8-bit path out of reach. +/// +/// MTL5 accepts every pairing at the `dot` level and swaps the operands onto +/// whichever the machine implements, because a dot product is symmetric. The +/// asymmetry is real one level down -- `simd::reduce_dot_widen` rejects +/// `(int8, uint8)` at compile time -- but that restriction is the kernel's, and +/// re-exposing it here would refuse a call the library can serve. +template +double mixed_dot_pair(const VA& a, const VB& b, AccKind kind) { + if (a.size() != b.size()) + throw std::invalid_argument("dot: vectors must have the same length"); + nogil guard; + return dispatch_acc(kind, [&]() -> double { + if constexpr (std::is_same_v) + return static_cast(mtl::dot<>(a, b)); + else + return static_cast(mtl::dot(a, b)); + }); +} + +/// The narrow integer element types (i8, i16, u8): `dot` only, and only with an +/// accumulator wider than the element. +/// +/// This is what phase 1 registered the containers for. MTL5 v5.11.0 routes +/// `dot` over 8- and 16-bit operands onto the hardware widening +/// multiply-accumulate -- vpmaddwd / vpdpbusd on x86, SMLAL / SDOT on NEON -- +/// and an int32 accumulator is what those instructions accumulate into. +/// +/// `norm` and `frobenius_norm` are deliberately absent, for the reason phase 1 +/// gave: two_norm takes sqrt of the accumulated sum, and on a wrapped -- and for +/// i16 negative -- sum that yields nan. An int32 accumulator would fix the sum +/// but the API has no way to say "accumulate in int32, deliver a real square +/// root", so the honest surface is not to offer it yet. +/// +/// `accumulator` is REQUIRED here, unlike every other dtype. Omitting it means +/// element precision, which for these widths is exactly the wrapping that phase +/// 1 refused to expose: six products of 100 give 96 for i8 where the answer is +/// 60000. Rather than silently redefine None for these dtypes alone, ask. +/// Register `dot` for one ordered pair of distinct narrow integer operands. +template +void register_mixed_narrow_int_pair(nb::module_& mx) { + const char* kPairDoc = + "Dot product over mixed-signedness 8-bit operands, accumulated in int32.\n\n" + "accumulator='i32' is required. u8 x i8 is VNNI's native pairing on x86\n" + "(unsigned activations against signed weights); ARM implements the\n" + "symmetric pairings first instead. Either order is accepted here -- a dot\n" + "product is symmetric, so MTL5 swaps the operands onto whichever form the\n" + "machine has rather than dropping to the generic loop.\n\n" + "The overflow contract is the same as the same-type form: products are\n" + "exact, the sum wraps, and 8-bit operands give roughly 131000 terms of\n" + "headroom in an int32."; + + auto acc_of = [](const std::optional& spec) { + if (!spec) + throw std::invalid_argument( + "mixed.dot on mixed 8-bit operands requires an explicit " + "accumulator: pass accumulator='i32'."); + return parse_acc(spec, type_suffix(), /*quire_ok=*/false, /*i32_ok=*/true); + }; + + mx.def("dot", [acc_of](const VectorView& a, const VectorView& b, + std::optional accumulator, + std::optional result) { + if (result) + throw std::invalid_argument( + "mixed.dot on narrow integer operands does not accept result=."); + return mixed_dot_pair(a.vec, b.vec, acc_of(accumulator)); + }, "a"_a, "b"_a, "accumulator"_a = nb::none(), "result"_a = nb::none(), kPairDoc); + + mx.def("dot", [acc_of](nb::ndarray, nb::c_contig, nb::device::cpu> a, + nb::ndarray, nb::c_contig, nb::device::cpu> b, + std::optional accumulator, + std::optional result) { + if (result) + throw std::invalid_argument( + "mixed.dot on narrow integer operands does not accept result=."); + const AccKind kind = acc_of(accumulator); + const std::size_t n = a.shape(0); + auto va = mtl::vec::dense_vector(n, const_cast(a.data())); + auto vb = mtl::vec::dense_vector(n, const_cast(b.data())); + return mixed_dot_pair(va, vb, kind); + }, nb::arg("a").noconvert(), nb::arg("b").noconvert(), + "accumulator"_a = nb::none(), "result"_a = nb::none(), kPairDoc); +} + +template +void register_mixed_narrow_int(nb::module_& mx) { + using VV = VectorView; + record_quire_support(); + + auto acc_of = [](const std::optional& spec) { + if (!spec) + throw std::invalid_argument( + std::string("mixed.dot on '") + type_suffix() + "' requires an " + "explicit accumulator: pass accumulator='i32'. The default is " + "element precision, and an 8- or 16-bit accumulator overflows " + "almost immediately -- six products of 100 wrap to 96 where the " + "exact answer is 60000. See mtl5.mixed.accumulators('" + + type_suffix() + "')."); + return parse_acc(spec, type_suffix(), /*quire_ok=*/false, /*i32_ok=*/true); + }; + + // The overflow contract, stated where a caller will meet it. This is not + // decoration: the headroom differs by four orders of magnitude across the + // operand widths, and a caller who does not know that gets silent wraparound. + const char* kDoc = + "Dot product accumulated in a precision wider than the element.\n\n" + "accumulator='i32' is required. On 8-bit operands this is the quad\n" + "multiply-accumulate (vpdpbusd / SDOT); on 16-bit it is the widening\n" + "multiply-accumulate (vpmaddwd / SMLAL).\n\n" + "OVERFLOW IS THE CONTRACT, not an error case. Products are always\n" + "exact -- two int16 cannot overflow an int32 product -- but the SUM\n" + "wraps, and how soon depends on the operand magnitude rather than the\n" + "vector length: at b bits of magnitude the headroom is about 2^(31-2b)\n" + "terms. Measured at full range: one i16 x i16 product uses 2^30 of the\n" + "int32 range, so TWO of them already overflow it, while i8 x i8 holds\n" + "131071 terms. That five order of magnitude gap is why the\n" + "quantized-inference instructions are 8-bit. Wrapping is two's\n" + "complement and therefore bit-identical across lane counts, backends\n" + "and thread partitions -- reproducible, but still wrapping.\n\n" + "Every signedness pairing is accepted. The hardware implements only\n" + "some of them (x86 does unsigned x signed first, ARM the symmetric\n" + "ones), and a dot product is symmetric, so MTL5 swaps the operands onto\n" + "whichever the machine has rather than refusing or falling back to the\n" + "generic loop."; + + mx.def("dot", [acc_of](const VV& a, const VV& b, + std::optional accumulator, + std::optional result) { + if (result) + throw std::invalid_argument( + std::string("mixed.dot on '") + type_suffix() + "' does not accept " + "result=: rounding an int32 sum back to an 8- or 16-bit element " + "would re-introduce the wrap the accumulator exists to avoid. The " + "result is delivered as a Python float, which holds every int32 " + "exactly."); + return mixed_dot(a.vec, b.vec, acc_of(accumulator), /*result_element=*/false); + }, "a"_a, "b"_a, "accumulator"_a = nb::none(), "result"_a = nb::none(), kDoc); + + mx.def("dot", [acc_of](nb::ndarray, nb::c_contig, nb::device::cpu> a, + nb::ndarray, nb::c_contig, nb::device::cpu> b, + std::optional accumulator, + std::optional result) { + if (a.shape(0) != b.shape(0)) + throw std::invalid_argument("dot: vectors must have the same length"); + if (result) + throw std::invalid_argument( + std::string("mixed.dot on '") + type_suffix() + "' does not accept " + "result=: see the docstring."); + const AccKind kind = acc_of(accumulator); + const std::size_t n = a.shape(0); + auto va = mtl::vec::dense_vector(n, const_cast(a.data())); + auto vb = mtl::vec::dense_vector(n, const_cast(b.data())); + return mixed_dot(va, vb, kind, /*result_element=*/false); + }, nb::arg("a").noconvert(), nb::arg("b").noconvert(), + "accumulator"_a = nb::none(), "result"_a = nb::none(), kDoc); +} + // =========================================================================== // convert() — element-wise re-quantization into a target number system // @@ -490,6 +667,16 @@ void register_mixed_precision(nb::module_& m) { register_mixed_native(mx); register_mixed_native(mx); + // Narrow integers: dot only, accumulator required. See + // register_mixed_narrow_int for why norm/frobenius_norm are absent. + register_mixed_narrow_int(mx); + register_mixed_narrow_int(mx); + register_mixed_narrow_int(mx); + // Mixed-signedness 8-bit pairings, both orders. u8 x i8 is what VNNI + // implements natively on x86 and what quantized inference is written in. + register_mixed_narrow_int_pair(mx); + register_mixed_narrow_int_pair(mx); + register_mixed_universal(mx); register_mixed_universal(mx); register_mixed_universal(mx); @@ -513,6 +700,12 @@ void register_mixed_precision(nb::module_& m) { register_mixed_universal(mx); mx.def("accumulators", [](const std::string& dtype) { + // The narrow integer types answer differently: an int32 accumulator is + // the one that maps to hardware, and the float accumulators are not + // offered because dot over int8 operands takes the generic loop + // and would quietly be the slow path dressed as a precision choice. + if (dtype == "i8" || dtype == "i16" || dtype == "u8") + return std::vector{"i32"}; std::vector v{"f32", "f64", "fma32", "fma64"}; // Answered from the registry the registration templates fill with // quire_for::ok -- the SAME trait dispatch_acc() consults -- so what @@ -554,7 +747,17 @@ void register_mixed_precision(nb::module_& m) { "fixpnt8", "fixpnt16", "lns16", "lns32", "cfloat32", "takum32", "dd_cascade", "td_cascade", "qd_cascade"}; - }, "Element dtypes accepted by convert() and the mixed-precision operations"); + }, "Element dtypes accepted by convert() and the mixed-precision operations.\n\n" + "This is the set convert() can TARGET, and callers rely on that -- the\n" + "test suite parametrizes over it and converts into every entry.\n\n" + "The narrow integer element types (i8, i16, u8) are deliberately NOT\n" + "here. They are real element types with containers and a mixed.dot\n" + "(accumulator='i32'), but convert() cannot target them: it re-quantizes\n" + "a float64 array, and rounding reals into an 8-bit integer is a\n" + "quantization scheme -- scale, zero point, rounding mode -- rather than\n" + "a cast. A naive version would silently clip everything outside\n" + "[-128, 127]. Build them with NumPy and pass them in;\n" + "mtl5.mixed.accumulators('i8') answers for them."); // ----- Dense mixed-precision iterative refinement ------------------------- mx.def("lu_iterative_refine", diff --git a/python/src/mtl5_types.hpp b/python/src/mtl5_types.hpp index a519bbd..e5b0043 100644 --- a/python/src/mtl5_types.hpp +++ b/python/src/mtl5_types.hpp @@ -228,20 +228,33 @@ void register_dense_ops(nb::module_& m); // Shared between the mixed-precision operations and the sparse factorizations // so that `accumulator=` means the same thing everywhere it appears. // =========================================================================== -enum class AccKind { Default, F32, F64, FMA32, FMA64, Quire }; +enum class AccKind { Default, F32, F64, FMA32, FMA64, Quire, I32 }; inline constexpr const char* kAccumulatorHelp = "valid accumulators: None (element precision), 'f32', 'f64', " - "'fma32', 'fma64'/'fma', 'quire'"; + "'fma32', 'fma64'/'fma', 'quire', 'i32' (8- and 16-bit integer operands)"; +/// `i32_ok` is true only for the narrow integer element types, where an int32 +/// accumulator is the whole point -- it is what vpmaddwd / vpdpbusd / SDOT +/// accumulate into. Offering it on a float or posit dtype would be meaningless, +/// so it is rejected by name rather than silently ignored. inline AccKind parse_acc(const std::optional& spec, - const char* dtype, bool quire_ok) { + const char* dtype, bool quire_ok, bool i32_ok = false) { if (!spec || *spec == "none" || *spec == "default") return AccKind::Default; const std::string& a = *spec; if (a == "f32" || a == "float32") return AccKind::F32; if (a == "f64" || a == "float64") return AccKind::F64; if (a == "fma32") return AccKind::FMA32; if (a == "fma" || a == "fma64") return AccKind::FMA64; + if (a == "i32" || a == "int32") { + if (!i32_ok) + throw std::invalid_argument( + std::string("accumulator='i32' is not available for dtype '") + dtype + + "': an int32 accumulator is for 8- and 16-bit INTEGER operands " + "(i8, i16, u8), which is what the hardware widening " + "multiply-accumulate takes. Use 'f64' to accumulate in double."); + return AccKind::I32; + } if (a == "quire") { if (!quire_ok) throw std::invalid_argument( diff --git a/tests/test_mixed_precision.py b/tests/test_mixed_precision.py index 406644c..155d5a8 100644 --- a/tests/test_mixed_precision.py +++ b/tests/test_mixed_precision.py @@ -50,6 +50,18 @@ def exact_sum_of_squares(a) -> float: ] +def _wrap32(x: int) -> int: + """Reduce an exact integer to what an int32 accumulator holds. + + Two's-complement wrapping is MTL5's documented contract for integer lanes, + and it is what makes an integer reduction bit-identical across lane counts, + backends and thread partitions. Tests assert against this rather than + against "no overflow", so they stay true in the regime these operand widths + actually run in. + """ + return ((x + 2**31) % 2**32) - 2**31 + + class TestConvert: def test_dtypes_listed(self): d = mtl5.dtypes() @@ -420,3 +432,118 @@ def test_refinement_through_ilu0_now_converges(self): ) assert info["converged"] assert np.linalg.norm(x - xt) / np.linalg.norm(xt) < 1e-10 + + +class TestIntegerAccumulator: + """`accumulator='i32'` over 8- and 16-bit integer operands (#88 phase 2). + + MTL5 v5.11.0 routes `dot` over narrow integer operands onto the + hardware widening multiply-accumulate — vpmaddwd / vpdpbusd on x86, SMLAL / + SDOT on NEON. An int32 accumulator is what those instructions accumulate + into, so it is the accumulator this path exists for. + """ + + NARROW = [(np.int8, "i8"), (np.int16, "i16"), (np.uint8, "u8")] + + @pytest.mark.parametrize("dt,suffix", NARROW, ids=[s for _, s in NARROW]) + def test_accumulators_offers_only_i32(self, dt, suffix): + assert mtl5.mixed.accumulators(suffix) == ["i32"] + + @pytest.mark.parametrize("dt,suffix", NARROW, ids=[s for _, s in NARROW]) + def test_matches_an_exact_int64_reference(self, dt, suffix): + rng = np.random.default_rng(0) + info = np.iinfo(dt) + a = rng.integers(info.min, info.max + 1, size=2000, dtype=dt) + b = rng.integers(info.min, info.max + 1, size=2000, dtype=dt) + got = mtl5.mixed.dot(a, b, accumulator="i32") + want = int(np.dot(a.astype(np.int64), b.astype(np.int64))) + assert int(got) == _wrap32(want) + + @pytest.mark.parametrize("dt,suffix", NARROW, ids=[s for _, s in NARROW]) + def test_container_and_ndarray_forms_agree(self, dt, suffix): + a = np.arange(1, 17, dtype=dt) + via_np = mtl5.mixed.dot(a, a, accumulator="i32") + via_vec = mtl5.mixed.dot(mtl5.vector(a), mtl5.vector(a), accumulator="i32") + assert via_np == via_vec + + @pytest.mark.parametrize("dt,suffix", NARROW, ids=[s for _, s in NARROW]) + def test_accumulator_is_required(self, dt, suffix): + """Omitting it means element precision, which is exactly the wrapping + phase 1 refused to expose. Ask rather than silently redefine None.""" + v = mtl5.vector(np.full(6, 100, dtype=dt)) + with pytest.raises(ValueError, match="requires an explicit accumulator"): + mtl5.mixed.dot(v, v) + + @pytest.mark.parametrize("dt,suffix", NARROW, ids=[s for _, s in NARROW]) + def test_result_element_is_refused(self, dt, suffix): + """Rounding the int32 sum back to an 8- or 16-bit element would + re-introduce the wrap the accumulator exists to avoid.""" + v = mtl5.vector(np.full(6, 100, dtype=dt)) + with pytest.raises(ValueError, match="does not accept result="): + mtl5.mixed.dot(v, v, accumulator="i32", result="element") + + def test_i32_is_refused_on_non_integer_dtypes(self): + v = mtl5.vector(np.ones(4)) + with pytest.raises(ValueError, match="not available for dtype"): + mtl5.mixed.dot(v, v, accumulator="i32") + + def test_the_wrap_is_twos_complement_not_undefined(self): + """Overflow is the contract, not an error case — and it is well defined. + + This is what the v5.11.0 pin buys: before it, the generic integer loops + were UB on overflow rather than the documented modular wrap, and for + these operand widths overflow is the normal regime. + + Full range i8: each -128 * -128 product is 16384, so an int32 holds + exactly 131071 terms and 131072 is the first that does not. + """ + for n in (131071, 131072, 131073): + a = np.full(n, -128, dtype=np.int8) + exact = n * 16384 + assert int(mtl5.mixed.dot(a, a, accumulator="i32")) == _wrap32(exact) + assert 131071 * 16384 <= 2**31 - 1, "131071 terms must still be exact" + assert 131072 * 16384 > 2**31 - 1, "...and 131072 must not" + + def test_i16_headroom_is_two_products(self): + """One full-range i16 product uses 2^30 of the int32 range, so two + already overflow it — five orders of magnitude less headroom than i8, + which is why the quantized-inference instructions are 8-bit.""" + one = np.full(1, -32768, dtype=np.int16) + assert int(mtl5.mixed.dot(one, one, accumulator="i32")) == 2**30 + two = np.full(2, -32768, dtype=np.int16) + assert int(mtl5.mixed.dot(two, two, accumulator="i32")) == _wrap32(2**31) + + +class TestIntegerAccumulatorMixedSignedness: + """`u8 x i8` is VNNI's native pairing on x86 — unsigned activations against + signed weights — and is what quantized inference is written in. + + MTL5 accepts every pairing at the `dot` level and swaps the operands onto + whichever form the machine implements, because a dot product is symmetric. + The kernel below it (`simd::reduce_dot_widen`) rejects `(int8, uint8)`, but + that restriction is the kernel's; re-exposing it here would refuse a call + the library can serve. + """ + + def _operands(self): + rng = np.random.default_rng(7) + u = rng.integers(0, 256, size=2048, dtype=np.uint8) + i = rng.integers(-128, 128, size=2048, dtype=np.int8) + return u, i, int(np.dot(u.astype(np.int64), i.astype(np.int64))) + + def test_u8_times_i8(self): + u, i, exact = self._operands() + assert int(mtl5.mixed.dot(u, i, accumulator="i32")) == _wrap32(exact) + + def test_i8_times_u8_is_accepted_too(self): + u, i, exact = self._operands() + assert int(mtl5.mixed.dot(i, u, accumulator="i32")) == _wrap32(exact) + + def test_the_two_orders_agree(self): + u, i, _ = self._operands() + assert mtl5.mixed.dot(u, i, accumulator="i32") == mtl5.mixed.dot(i, u, accumulator="i32") + + def test_container_form_works_for_pairs(self): + u, i, exact = self._operands() + got = mtl5.mixed.dot(mtl5.vector(u), mtl5.vector(i), accumulator="i32") + assert int(got) == _wrap32(exact)