From de5acfd8c671fa404c9f1e030d07e52634080971 Mon Sep 17 00:00:00 2001 From: Ettore Di Giacinto Date: Tue, 18 Aug 2026 09:50:33 +0000 Subject: [PATCH 1/3] spec(VT-MATMUL-FP8-BLOCK-REF): the block-scaled FP8 GEMM applies its scales in the MAINLOOP, which an epilogue alpha cannot express Milestone M2 of #1189, specified before any implementation exists. The op is `vt::MatmulFp8BlockScaled`, the CPU reference arm of the 128x128 block-scaled FP8 GEMM, and it is the numerical oracle every later kernel is measured against. M1 (`ad5f175e7`) established that upstream's Triton source is not what executes and that the two arms disagree in polarity, so the same question was asked here before any arithmetic was mirrored. The answer runs through six files and it lands differently: the executing kernel is CUTLASS, and it AGREES with the reference. `vllm/model_executor/kernels/linear/__init__.py:355-377` ranks CUTLASS third; DeepGEMM is auto-disabled for `qwen3_5_text` on family 120 (`vllm/utils/deep_gemm.py:27-46`); `scaled_mm_entry.cu:220-226` routes every `sm >= 120` device to the sm120 entry; `scaled_mm_helper.hpp:39-55` takes the blockwise branch and checks both shapes with `ceil_div`; and `scaled_mm_blockwise_sm120_fp8_dispatch.cuh:218-235` hands both scale pointers to the MAINLOOP arguments over a `float` accumulator. CUTLASS 4.5.0's mainloop line is `accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i)` (`sm120_mma_tma_blockwise_scaling.hpp:714-717`), which is `native_w8a8_block_matmul`'s `c += matmul(a, b.t()) * s` (`tests/kernels/quant_utils.py:145-151`) with the two scale multiplies associated left to right instead of as a product. Upstream's own test compares the two at `rel_diff < 0.001` (`test_block_fp8.py:194-200`), which is what admits that one-ULP difference and what no polarity error would survive. The constraint the row exists for: our per-tensor FP8 path folds one scalar `alpha` in the epilogue, and an epilogue has exactly one degree of freedom per output element while the block scheme has `cdiv(K, block_k)` of them. An epilogue-only application therefore cannot express a per-K-block scale AT ALL. That is a correctness constraint rather than an optimisation choice, and it is why this is a new op rather than a parameter of `kMatmulFp8Cutlass`. G4 is the instrument: two K-blocks with different scales, then the same operands with those two scales SWAPPED, which leaves every per-tensor summary of the scales identical and changes the correct answer. Ragged edges are legal and gated. The production wrapper asserts `triton.cdiv(N, block_n) == Bs.shape[0]` and `triton.cdiv(K, block_k) == Bs.shape[1]` (`fp8_utils.py:935-936`), so a short final block must work. `N=576` and `K=3884` from upstream's grid are the shapes that expose an integer-division bug and both are in the ported grid, separately and together. The upstream grid is adapted and the adaptation is stated: 125 combinations whose largest is about 9.3e11 multiply-accumulates is a GPU grid, and a naive CPU reference nest run twice cannot execute it. Every axis VALUE is preserved, the pairing is not, and the parameters, dtypes, tolerance and failure criterion are preserved exactly. M2 lands unreached, which is named here rather than discovered later. No production entry point will dispatch `vt::MatmulFp8BlockScaled` at its merge commit: `include/vllm.h` does not expose it, no loader builds an `Fp8BlockWeight`, and `ModelRegistry::Forward` has no block-FP8 linear method. Milestone M4 owns the wiring and needs M3 first. The spec lists it under `## Owed`, per the staged-slice exception in `.agents/reachability.md`. No `.agents/issue-index.md` row is appended: M1 already added the #1189 row and that file is append-only, so re-stating the issue would duplicate rather than merge. FOLLOWING_AGENTS_PROTOCOL Following-Agents-Protocol: true AI-Assisted: true Assisted-by: AGENT:claude-opus-5-1m [claude-code] --- .agents/specs/vt-matmul-fp8-block-ref.md | 296 +++++++++++++++++++++++ 1 file changed, 296 insertions(+) create mode 100644 .agents/specs/vt-matmul-fp8-block-ref.md diff --git a/.agents/specs/vt-matmul-fp8-block-ref.md b/.agents/specs/vt-matmul-fp8-block-ref.md new file mode 100644 index 000000000..1d5195e6c --- /dev/null +++ b/.agents/specs/vt-matmul-fp8-block-ref.md @@ -0,0 +1,296 @@ +# VT-MATMUL-FP8-BLOCK-REF — the block-scaled FP8 GEMM, CPU reference arm + +Issue: [#1189](https://github.com/mudler/vllm.cpp/issues/1189), milestone **M2**. +Row: `VT-MATMUL-FP8-BLOCK-REF`. +Pinned oracle: vLLM `5559679229bc961848b121ccdeaa8fa5d79bec98` +(`.agents/upstream-sync.md`), asserted as the HEAD of the local checkout before +any `file:line` below was read. + +## Scope + +Add one `vt` op with a CPU kernel only: + +```c++ +void MatmulFp8BlockScaled(Queue& q, Tensor& out, + const Tensor& a_fp8, const Tensor& a_scale, + const Tensor& b_fp8, const Tensor& b_scale, + int block_n, int block_k); +``` + +`a_fp8` is `[M,K]` i8 carrying raw fp8-e4m3fn bytes, as `vt::QuantFp8Group` +emits. `a_scale` is `[M, cdiv(K, block_k)]` f32, as `vt::QuantFp8Group` emits. +`b_fp8` is `[N,K]` i8, the on-disk weight. `b_scale` is +`[cdiv(N, block_n), cdiv(K, block_k)]` f32, the on-disk `weight_scale_inv`. +`out` is `[M,N]` f32 or bf16. + +This is the numerical oracle every later block-FP8 kernel is measured against. +It is milestone M2 of #1189 and it deliberately stops there. + +**Out of scope, each owned by another milestone of #1189**: `Fp8BlockWeight`, +the loader rung and the config reader (M3); `layers::Fp8BlockLinearMethod` and +the Qwen3.5 wiring (M4); the mainloop-scaled CUTLASS kernel for `sm_121a` (M5); +merged `gate_up` and QKV (M6). No CUDA arm lands here: M5 owns it, and a +device-free host cannot gate one. + +## Which implementation actually runs, and why it does not change the answer + +M1 established that upstream's Triton source is not what executes for the +activation quant, and that the two arms disagree in polarity. The same question +has to be asked here before any arithmetic is mirrored, and the executing chain +was read end to end rather than inferred from the Python. + +The dispatch decision, in order: + +| Step | Where | +|---|---| +| CUDA block-FP8 kernel priority list: FlashInfer, DeepGEMM, **CUTLASS**, Marlin, Triton, Humming | `vllm/model_executor/kernels/linear/__init__.py:355-377` | +| DeepGEMM is auto-disabled for `qwen3_5_text` on device-capability family 120 | `vllm/utils/deep_gemm.py:27-46` | +| the selected kernel's apply: `ops.cutlass_scaled_mm(A, B.T, scale_a=As, scale_b=Bs.T)` | `vllm/model_executor/kernels/linear/scaled_mm/cutlass.py:312-326` | +| `cutlass_scaled_mm` routes every `sm >= 120` device to the sm120 entry | `csrc/libtorch_stable/quantization/w8a8/cutlass/scaled_mm_entry.cu:220-226` | +| the sm120 entry hands `cutlass_scaled_mm_blockwise_sm120_fp8` to the shared dispatcher | `csrc/libtorch_stable/quantization/w8a8/cutlass/scaled_mm_c3x_sm120.cu:13-20` | +| the dispatcher takes the blockwise branch whenever a scale is not per-tensor or per-row, requires both scales **f32** and 2-D, checks the shapes with **`ceil_div`**, and refuses a bias | `csrc/libtorch_stable/quantization/w8a8/cutlass/c3x/scaled_mm_helper.hpp:15-18,39-55` | +| the blockwise kernel's accumulator and scale element type are both `float`, and the two scale pointers are **mainloop** arguments | `csrc/libtorch_stable/quantization/w8a8/cutlass/c3x/scaled_mm_blockwise_sm120_fp8_dispatch.cuh:56-58,218-235` | +| CUTLASS 4.5.0, the mainloop line itself: `accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i)` | `include/cutlass/gemm/collective/sm120_mma_tma_blockwise_scaling.hpp:714-717` | + +So the executing kernel is CUTLASS, and unlike M1 it **agrees** with the +reference in both placement and polarity: a per-K-block partial product is +formed in an f32 accumulator and **multiplied** by `a_s * b_s` before it is +added to the running sum. Mirroring `native_w8a8_block_matmul` is therefore +mirroring the executing arm, and that is a measurement rather than an +assumption: upstream compares the CUTLASS arm against exactly this reference at +`rel_diff < 0.001` (`tests/kernels/quantization/test_block_fp8.py:194-200`), a +gate no polarity or placement error survives. + +One difference is real and is not a divergence. The reference forms the scale +product first, `s = As_tiles[i] * Bs[j][i]` then `partial * s` +(`tests/kernels/quant_utils.py:150-151`); CUTLASS associates left to right, +`(partial * a_s) * b_s` (`sm120_mma_tma_blockwise_scaling.hpp:715`). The two +differ by at most one f32 ULP per K-block and upstream's own 1e-3 relative gate +is what admits it. We mirror the reference's association, because this op **is** +the reference port and M5's kernel will be measured against it. + +## The constraint that decides correctness + +**The scales apply in the mainloop, once per K-block, into an f32 accumulator — +not in the epilogue.** Written out (`fp8_utils.py:826-836`, which is the Triton +arm spelling the same structure the CUTLASS mainloop implements): + +```text +accumulator = zeros(f32) +for k_block: + a_s = a_scale[m, k_block] per token, per K-group + b_s = b_scale[n_block, k_block] per (N-block, K-block) + accumulator += dot(a_tile, b_tile) * a_s * b_s +``` + +Our existing per-tensor FP8 path folds one scalar `alpha` into the epilogue +(`vt::MatmulFp8Cutlass`, `src/vt/cpu/cpu_ops.cpp` `MatmulFp8CutlassKernel`: +`alpha * acc` after the whole K reduction). **An epilogue-only application +cannot express a per-K-block scale at all**: it has exactly one degree of +freedom per output element, and the block scheme has `cdiv(K, block_k)` of them. +This is a correctness constraint, not an optimisation choice, and it is the +reason M2 is a new op rather than a parameter of `kMatmulFp8Cutlass`. + +G4 is the instrument for it, and it is constructed so that no epilogue-folded +alpha can pass: it builds a case whose two K-blocks carry different scales and +whose correct output is not `alpha * (plain fp8 GEMM)` for any `alpha`. + +## Upstream anchors + +| What | Where | +|---|---| +| the reference this op ports, whole | `tests/kernels/quant_utils.py:91-154` | +| the fp32 compute type and the `.to(output_dtype)` store | `quant_utils.py:98,109-110,153` | +| `ceil` tiling of both dimensions | `quant_utils.py:123-124` | +| the `Bs` shape assertions, `n_tiles`/`k_tiles` | `quant_utils.py:125-126` | +| the `As` last-dimension `ceil` assertion | `quant_utils.py:115-116` | +| the mainloop accumulate, scales multiplied per K-block | `quant_utils.py:145-151` | +| the production wrapper's shape contract, `cdiv` on N and K | `vllm/model_executor/layers/quantization/utils/fp8_utils.py:928-936` | +| the Triton mainloop spelling the same structure | `fp8_utils.py:826-836` | +| the ported test, its grid, its tolerance | `tests/kernels/quantization/test_block_fp8.py:48-55,123-153` | +| the ragged-shape case upstream calls out by name (DSV3 `kv_a_proj_with_mqa`) | `test_block_fp8.py:156-200` | +| the activation quant that produces `a_fp8`/`a_scale` | `.agents/specs/vt-quant-fp8-group.md`, landed `ad5f175e7` | + +## Design + +### Numerics + +```text +for each output element (m, n): + acc = 0 f32 + for kt in [0, cdiv(K, block_k)): + part = 0 f32, a SEPARATE accumulator + for k in the k-tile (ragged final tile is short): + part += f32(a_fp8[m,k]) * f32(b_fp8[n,k]) + acc += part * (a_scale[m, kt] * b_scale[n / block_n, kt]) + out[m, n] = acc stored to out's dtype +``` + +`part` is a separate register that is scaled and then folded into `acc`. That is +the whole point: collapsing it into `acc` is what an epilogue alpha does, and it +is unrepresentable. + +`n / block_n` is integer division of the **output column index**, mirroring +`offs_bsn = offs_bn // group_n` (`fp8_utils.py:823`) and the reference's tiling +at `quant_utils.py:131-143`. `block_n` and `block_k` are validated positive +before either divides anything, because `x % 0` and `x / 0` are undefined +behaviour and a zero must refuse rather than trap. + +### Ragged edges + +`cdiv`, not floor, on both dimensions. Upstream's production wrapper asserts +`triton.cdiv(N, block_n) == Bs.shape[0]` and +`triton.cdiv(K, block_k) == Bs.shape[1]` (`fp8_utils.py:935-936`), so a final +short block is legal and must work. The final K-tile runs to `K`, not to +`(kt+1)*block_k`; the final N-block is short and its scale row still exists. +`N=576` (`4*128 + 64`) and `K=3884` (`30*128 + 44`) are the shapes from +upstream's own grid that expose an integer-division bug, and both are in the +ported grid, separately and together. + +### Memory format + +`out` is f32 or bf16 and nothing here is f32 by default. The accumulator is f32 +because upstream's `ElementAccumulator` is `float` +(`scaled_mm_blockwise_sm120_fp8_dispatch.cuh:56`) and the reference's +`compute_type` is `torch.float32` (`quant_utils.py:98`); the **store** rounds to +whatever `out` carries, exactly as `.to(output_dtype)` does +(`quant_utils.py:153`). Both scale tensors are f32 and the op refuses any other +dtype, because upstream refuses it too +(`scaled_mm_helper.hpp:15-18`). This op therefore widens nothing: it is the same +polarity `.agents/porting.md` requires. + +### Structure + +`vt::MatmulNvfp4Fp4` (`include/vt/ops.h:1477`) is the closest existing +signature — a block-scaled GEMM taking two packed operands and two scale +streams — and `MatmulFp8BlockScaled` follows its shape, with the scalar `alpha` +replaced by the block geometry it cannot express. `OpId::kMatmulFp8BlockScaled` +is appended before `kCount`, the additive convention documented at +`include/vt/ops.h:363-368`, so no existing op's id shifts. + +The kernel is a **correctness reference, not a performance path**, in the same +sense as `MatmulFp8CutlassKernel` beside it: a naive nest that makes the +block-FP8 seam resolvable on a CPU queue so M3 and M4 can be gated without a +GPU. It makes no speed claim. + +## Risks + +| Risk | Control | +|---|---| +| a later reader folds the scales into the epilogue, "simplifying" the inner accumulator away | the prose above and in the header, plus G4, which is unpassable by any single-alpha form | +| a ragged N or K uses floor tiling and silently drops or misindexes a block | `cdiv` everywhere; G2 runs `N=576`, `K=3884` separately and together, and G5 refuses a floor-sized scale tensor by name | +| the `n / block_n` index is computed from the tile index rather than the column, which agrees for round N and diverges for ragged N | G2's `N=576` case, where the two disagree in the final block | +| an all-zero output passes every value comparison | G2, G3 and G4 each carry a `nonzero == numel` vacuity guard | +| a scale dtype silently widens or narrows | G5 refuses a non-f32 scale by name | +| the test's reference is the implementation rewritten, so agreement proves nothing | the G2/G3 reference accumulates in `double` with a different loop nest and an exhaustively derived fp8 decode table, and its bound is an f32 forward-error bound rather than a fudge factor | +| `block_n`/`block_k` of zero divides by zero | validated positive before use; G5 asserts the refusal | +| nothing dispatches the op, so it is dead on arrival | acknowledged and named under `## Owed`: this is the staged-slice exception of `.agents/reachability.md`, and M4 owns the wiring | + +## Tests + +`tests/vt/test_ops_matmul_fp8_block_cpu.cpp`, registered in +`tests/CMakeLists.txt`. + +- **G1** the registration itself: `OpRegistered(kMatmulFp8BlockScaled, kCPU)`, + and `OpName` is not `unknown`. A refusal case cannot stand in for this — M1 + measured a deleted registration passing every refusal test it had. +- **G2** the ported upstream case, `test_block_fp8.py:123-153`, against an + independently written `double` reference, carrying upstream's own + `rel_diff < 0.001` criterion verbatim **and** a tighter per-element f32 + forward-error bound, because our arm is the reference rather than a kernel + measured against it. Vacuity guard. +- **G3** the ragged grid, run as part of G2's table and called out here because + it is the reason the grid is what it is: `N=576`, `K=3884`, and the pair + together, plus upstream's dedicated DSV3 case `M=32, N=576, K=7168` + (`test_block_fp8.py:156-200`). +- **G4** **the mainloop constraint, made unpassable by an epilogue**. Two + K-blocks with deliberately different scales and a hand-computed expected + output; then the same operands with the two K-block scales *swapped*, which + leaves every per-tensor summary of the scales identical and changes the + correct answer. A kernel that folds one alpha produces the same output for + both and fails. Vacuity guard. +- **G5** the refusals, each by name: a non-f32 scale, a `b_scale` sized by + floor instead of `cdiv` on N and on K, an `a_scale` with the wrong number of + groups or rows, a zero or negative `block_n`/`block_k`, a rank mismatch, a + non-contiguous operand, a device mismatch, a non-i8 packed operand, and an + `out` that is neither f32 nor bf16. +- **G6** the M1 seam: `vt::QuantFp8Group` feeding `vt::MatmulFp8BlockScaled` on + a CPU queue, end to end, which is the pair a block-FP8 linear method will run. + This is a *composition* test, not a reachability claim; nothing in production + calls it yet and `## Owed` says so. + +### The adaptation of the upstream grid, and why it is unavoidable + +Upstream's grid is `itertools.product(M, N, K, ...)` with `M=[1,7,8,83,4096]`, +`N=[128,512,576,7168,13824]`, `K=[256,3884,4096,13824,16384]` +(`test_block_fp8.py:48-55`) — 125 combinations, the largest of which is +`4096x13824x16384`, about 9.3e11 multiply-accumulates. That is a GPU grid. A +naive CPU reference nest cannot run it, and running it twice (op and reference) +is worse. + +The adaptation preserves the **axes**, not the product: every value of every +axis appears at least once, the ragged values appear separately and together, +and the total is under 1e9 MACs. Which values are paired with which is the only +thing dropped, and the parameters, dtypes, tolerance and failure criterion are +preserved exactly. The grid is written out in the test with this reasoning +beside it. + +## Gates + +| Gate | Command | +|---|---| +| focused | `ctest -R test_ops_matmul_fp8_block_cpu --output-on-failure` | +| op provider totality | `ctest -R test_op_provider` | +| the M1 sibling, unchanged | `ctest -R test_ops_quant_fp8_group_cpu` | +| the per-tensor sibling, unchanged | `ctest -R test_ops_fp8_cpu` | +| record | `scripts/agent-preflight.sh --fail-on-skip` | + +No GPU lease is taken and none is needed: the op has a CPU arm only. + +## Owed + +- **Nothing reaches this op.** `vt::MatmulFp8BlockScaled` is dispatched by no + production entry point at this merge commit: `include/vllm.h` does not expose + it, no loader builds an `Fp8BlockWeight`, and `ModelRegistry::Forward` has no + block-FP8 linear method to call it from. The wiring is owned by #1189 + milestone M4 (`layers::Fp8BlockLinearMethod` and the Qwen3.5 dense forward), + which needs M3 first. This is the staged-slice exception of + `.agents/reachability.md`, named here, in the commit body, and in the pull + request body. +- **The CUDA arm.** There is none. Milestone M5 owns the mainloop-scaled CUTLASS + kernel for `sm_121a`, needs a GPU, and will be measured against this op. +- **The `swap_ab` path.** Upstream's sm120 blockwise dispatch swaps the operands + for some shapes (`scaled_mm_blockwise_sm120_fp8_dispatch.cuh:221-235`), which + changes the scale-pointer assignment but not the mathematical result. It is a + kernel-scheduling concern with no reference-arm meaning. Owed by M5. +- **The column-major and TMA-aligned activation-scale layouts** + (`fp8_utils.py:610-628`), which the CUTLASS kernel reads and this reference + does not. Already owed by M5 in `.agents/specs/vt-quant-fp8-group.md`; repeated + here because M2 is the first op that could have consumed them. +- **`bias`.** Upstream refuses a bias on the blockwise path outright + (`scaled_mm_helper.hpp:55`). This op takes none, which mirrors the refusal + rather than deferring a feature. + +## Stop conditions + +Stop and report `NEEDS_DECISION` if any of the following holds. + +- The pinned oracle's checkout is not at + `5559679229bc961848b121ccdeaa8fa5d79bec98`. Every anchor above was read at + that revision, asserted before the first read. +- G4 cannot be made to fail an epilogue-folded implementation. If a single alpha + can reproduce the mainloop result on the constructed case, the case is wrong + and the constraint is untested, and widening the case is not the fix. +- The op cannot express upstream's contract without a parameter that no caller + in this tree passes. + +Stop and report `NEEDS_CONTEXT` if the work requires a GPU lease. The row is +scoped so that it does not. + +## Evidence + +Recorded below at completion: the RED before the implementation existed, the +GREEN after, the per-block gate counts, and the mutation table. Every mutation +prints `compile_rc` **and** `git diff --stat`, because a mutation that fails to +build and a mutation that never applied both read as a passing test, and M1 hit +each of those once. From 4e9e5d6513c7156bc7006d2ebb6a867e0bbd053d Mon Sep 17 00:00:00 2001 From: Ettore Di Giacinto Date: Tue, 18 Aug 2026 10:09:45 +0000 Subject: [PATCH 2/3] feat(VT-MATMUL-FP8-BLOCK-REF): add vt::MatmulFp8BlockScaled, the block-scaled fp8 GEMM whose scales cannot live in an epilogue MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Milestone **M2** of #1189: the CPU reference arm of the 128x128 block-scaled FP8 GEMM, and the numerical oracle every later block-fp8 kernel is measured against. It takes `a_fp8 [M,K]` i8 with `a_scale [M, cdiv(K, block_k)]` f32 — exactly what `vt::QuantFp8Group` emits, landed as M1 in `ad5f175e7` — plus `b_fp8 [N,K]` i8 with `b_scale [cdiv(N, block_n), cdiv(K, block_k)]` f32, and writes `out [M,N]` f32 or bf16. Spec: [`.agents/specs/vt-matmul-fp8-block-ref.md`](.agents/specs/vt-matmul-fp8-block-ref.md), committed before this change. Pinned oracle: vLLM `5559679229bc961848b121ccdeaa8fa5d79bec98`, asserted as the local checkout's HEAD before any anchor below was read. ## The constraint that decides correctness The scales apply in the GEMM MAINLOOP, once per K-block, into an f32 accumulator — not in the epilogue: ``` acc = 0 for kt in [0, cdiv(K, block_k)): part = 0 a SEPARATE register for k in the k-tile: part += f8(a[m,k]) * f8(b[n,k]) acc += part * ( a_scale[m,kt] * b_scale[n / block_n, kt] ) ``` `vt::MatmulFp8Cutlass` folds ONE scalar alpha after the whole K reduction. An epilogue has exactly one degree of freedom per output element; this scheme has `cdiv(K, block_k)` of them. An epilogue-only application therefore cannot express a per-K-block scale AT ALL, which is why this is a new op rather than a parameter of the existing one, and why `part` stays a separate register. ## Which implementation actually runs M1 found that upstream's Triton source is not what executes and that the two arms disagree in polarity, so the same question was asked here rather than assumed. The chain: `vllm/model_executor/kernels/linear/__init__.py:355-377` ranks CUTLASS third behind FlashInfer and DeepGEMM; DeepGEMM is auto-disabled for `qwen3_5_text` on device-capability family 120 (`vllm/utils/deep_gemm.py:27-46`) and Marlin is excluded at `cc >= 89`; `scaled_mm_entry.cu:220-226` routes every `sm >= 120` device to the sm120 entry; `scaled_mm_helper.hpp:39-55` takes the blockwise branch, requires both scales f32 and 2-D, checks the shapes with `ceil_div` and refuses a bias; and `scaled_mm_blockwise_sm120_fp8_dispatch.cuh:218-235` hands both scale pointers to the MAINLOOP arguments over an `ElementAccumulator = float`. CUTLASS 4.5.0's line is `accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i)` (`sm120_mma_tma_blockwise_scaling.hpp:714-717`). So the executing kernel is CUTLASS and, unlike M1, it AGREES with the reference: same placement, same polarity. The one difference is association — CUTLASS multiplies left to right where `native_w8a8_block_matmul` forms the scale product first (`tests/kernels/quant_utils.py:150-151`) — worth at most one f32 ULP per K-block, which is what upstream's own `rel_diff < 0.001` between the two admits (`test_block_fp8.py:194-200`). We mirror the reference's association, because this op IS the reference port and M5's kernel will be measured against it. ## Ragged edges `cdiv`, never floor, on both dimensions: upstream asserts `triton.cdiv(N, block_n) == Bs.shape[0]` and `triton.cdiv(K, block_k) == Bs.shape[1]` (`fp8_utils.py:935-936`), so a short final block is legal and must work. `n / block_n` indexes the b_scale row by OUTPUT COLUMN (`fp8_utils.py:823`), which agrees with a tile counter for round N and does not for ragged N. `N=576` and `K=3884` from upstream's grid are in the ported grid separately and together. ## Tests `tests/vt/test_ops_matmul_fp8_block_cpu.cpp`, registered in `tests/CMakeLists.txt`. RED first: before the op existed the focused build failed with 29 errors, 23 `'MatmulFp8BlockScaled' is not a member of 'vt'` and 6 `'kMatmulFp8BlockScaled' is not a member of 'vt::OpId'`. GREEN after: 6 cases, 80 assertions, 0 failed, in 7 s. Per block: G1 2, G2 22, G3 13, G4 14, G5 25, G6 4 — which sums to 80, so no bucket is silently empty. G4 is the instrument for the constraint above and is built so that no epilogue-folded alpha can pass it. Two K-blocks whose partials are 128 and 256 and whose scales are 0.25 and 0.5 give 160; SWAPPING the two K-block scales leaves every per-tensor summary of the scale tensor identical — same set, sum, product and max — and gives 128. Every value is exact in f32 and bf16, so these are equalities. A second subcase runs N=129, whose second N-block is one column wide, and checks each output against a hand-computed value. The upstream grid is adapted and the adaptation is written into the test. Upstream runs `itertools.product` over 125 combinations whose largest is `4096x13824x16384`, about 9.3e11 multiply-accumulates. That is a GPU grid and a naive CPU reference nest, run twice, cannot execute it. Every axis VALUE appears at least once and the ragged values appear separately and together; only the pairing is dropped, and the parameters, dtypes, tolerance and failure criterion are preserved exactly. The reference is independently written: upstream's loop nest, k-tile outer and n-block inner, in `double`, with an fp8 decode derived from the e4m3fn field layout rather than from any codec in `src/`. ## Nothing reaches this op yet `vt::MatmulFp8BlockScaled` is dispatched by no production entry point at this merge commit: `include/vllm.h` does not expose it, no loader builds an `Fp8BlockWeight`, and `ModelRegistry::Forward` has no block-FP8 linear method to call it from. The wiring is owned by #1189 milestone M4 (`layers::Fp8BlockLinearMethod` and the Qwen3.5 dense forward), which needs M3 first. This is the staged-slice exception of `.agents/reachability.md`, named here, in the pull request body, and under `## Owed` in the spec. G6 composes M1's quantizer with this GEMM end to end on a CPU queue, which is the pair a linear method will run; it is a composition test and not a reachability claim. No CUDA arm lands here. Milestone M5 owns the mainloop-scaled CUTLASS kernel for `sm_121a`, needs a GPU, and this row took no lease by design. FOLLOWING_AGENTS_PROTOCOL Following-Agents-Protocol: true AI-Assisted: true Assisted-by: AGENT:claude-opus-5-1m [claude-code] --- include/vt/ops.h | 78 +++ src/vt/cpu/cpu_ops.cpp | 76 +++ src/vt/op_provider.cpp | 2 + src/vt/ops.cpp | 44 ++ tests/CMakeLists.txt | 5 + tests/vt/test_ops_matmul_fp8_block_cpu.cpp | 684 +++++++++++++++++++++ 6 files changed, 889 insertions(+) create mode 100644 tests/vt/test_ops_matmul_fp8_block_cpu.cpp diff --git a/include/vt/ops.h b/include/vt/ops.h index 7905f4da4..a4142f6af 100644 --- a/include/vt/ops.h +++ b/include/vt/ops.h @@ -433,6 +433,14 @@ enum class OpId : uint8_t { // See vt::QuantFp8Group below for the contract. // Appended before kCount so no existing op's id shifts. kQuantFp8Group, + // --- Block-wise FP8 (VT-MATMUL-FP8-BLOCK-REF, #1189 milestone M2). The + // 128x128 block-scaled fp8 GEMM that consumes kQuantFp8Group's output. It is + // NOT a parameter of kMatmulFp8Cutlass and cannot be: that op folds ONE + // scalar alpha into the epilogue, which has exactly one degree of freedom per + // output element, while this scheme has cdiv(K, block_k) of them and applies + // them in the MAINLOOP. See vt::MatmulFp8BlockScaled below for the contract. + // Appended before kCount so no existing op's id shifts. + kMatmulFp8BlockScaled, kCount }; @@ -936,6 +944,13 @@ using QuantFp8StaticFn = void (*)(Queue&, Tensor&, const Tensor&, float); // and the f32 [M, K/group_size] per-group scale. using QuantFp8GroupFn = void (*)(Queue&, Tensor& /*out_fp8*/, Tensor& /*out_scale*/, const Tensor& /*x*/, int /*group_size*/); +// Two scale streams and the block geometry, where the per-tensor fp8 GEMMs above +// carry one scalar alpha. The geometry is not a convenience: it is what selects +// which scale pair each K-block multiplies by. +using MatmulFp8BlockScaledFn = void (*)(Queue&, Tensor& /*out*/, const Tensor& /*a_fp8*/, + const Tensor& /*a_scale*/, const Tensor& /*b_fp8*/, + const Tensor& /*b_scale*/, int /*block_n*/, + int /*block_k*/); using RmsNormQuantFp8Fn = void (*)(Queue&, Tensor& /*out_fp8*/, Tensor* /*out_bf16*/, const Tensor& /*x*/, const Tensor& /*weight*/, const RmsNormArgs&, Tensor* /*residual*/, float /*input_scale*/); @@ -1582,6 +1597,69 @@ void QuantFp8Static(Queue& q, Tensor& out_fp8, const Tensor& x, float input_scal void QuantFp8Group(Queue& q, Tensor& out_fp8, Tensor& out_scale, const Tensor& x, int group_size); +// MatmulFp8BlockScaled (VT-MATMUL-FP8-BLOCK-REF, #1189 M2, +// .agents/specs/vt-matmul-fp8-block-ref.md) — the 128x128 block-scaled fp8 GEMM, +// mirroring native_w8a8_block_matmul (tests/kernels/quant_utils.py:91-154): +// +// for each (m, n): +// acc = 0 f32 +// for kt in [0, cdiv(K, block_k)): +// part = 0 f32, a SEPARATE register +// for k in the k-tile: part += f8(a[m,k]) * f8(b[n,k]) +// acc += part * ( a_scale[m, kt] * b_scale[n / block_n, kt] ) +// out[m, n] = acc stored to out's dtype +// +// THE SCALES APPLY IN THE MAINLOOP, ONCE PER K-BLOCK, INTO AN F32 ACCUMULATOR — +// NOT IN THE EPILOGUE, and that is a correctness constraint rather than an +// optimisation choice. MatmulFp8Cutlass above folds one scalar alpha after the +// whole K reduction. An epilogue has exactly ONE degree of freedom per output +// element; this scheme has cdiv(K, block_k) of them. An epilogue-only +// application therefore cannot express a per-K-block scale AT ALL, which is why +// this is a separate op. DO NOT "simplify" `part` away into `acc`: that IS the +// epilogue form, and tests/vt/test_ops_matmul_fp8_block_cpu.cpp G4 is built so +// that no single-alpha implementation can pass it. +// +// WHICH UPSTREAM ARM THIS MIRRORS. The Triton kernel at fp8_utils.py:826-836 is +// not what executes on the target architecture; CUTLASS is +// (vllm/model_executor/kernels/linear/__init__.py:355-377 ranks it third, +// DeepGEMM is auto-disabled for qwen3_5_text on family 120 at +// vllm/utils/deep_gemm.py:27-46, and Marlin is excluded at cc >= 89). Unlike the +// QuantFp8Group case above, the two AGREE: csrc/.../c3x/ +// scaled_mm_blockwise_sm120_fp8_dispatch.cuh:56-58,218-235 hands both scale +// pointers to the MAINLOOP arguments over an `ElementAccumulator = float`, and +// cutlass 4.5.0's sm120_mma_tma_blockwise_scaling.hpp:714-717 is literally +// `accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i)`. +// CUTLASS associates the two scale multiplies left to right where the reference +// forms their product first (quant_utils.py:150-151); the difference is at most +// one f32 ULP per K-block and upstream's own gate admits it, comparing the two +// at rel_diff < 0.001 (test_block_fp8.py:194-200). We mirror the reference's +// association, because this op IS the reference port and the CUDA kernel that +// #1189 milestone M5 lands will be measured against it. +// +// SHAPES, with CEIL on every tiling, so a ragged final block is legal and must +// work (upstream asserts exactly this at fp8_utils.py:935-936): +// a_fp8 [M,K] i8, raw fp8-e4m3fn bytes +// a_scale [M, cdiv(K, block_k)] F32 +// b_fp8 [N,K] i8, raw fp8-e4m3fn bytes +// b_scale [cdiv(N, block_n), cdiv(K, block_k)] F32 +// out [M,N] f32 or bf16 +// The scales are f32 because upstream refuses any other dtype on this path +// (csrc/.../c3x/scaled_mm_helper.hpp:15-18) and the accumulator is f32. +// a_scale's K axis is a CEIL too, so this op accepts a K that vt::QuantFp8Group +// would refuse; that asymmetry is upstream's own (fp8_utils.py:930 uses cdiv +// where fp8_utils.py:596-599 demands divisibility). +// +// No bias: upstream refuses one outright on the blockwise path +// (scaled_mm_helper.hpp:54), so this mirrors a refusal rather than deferring a +// feature. +// +// CPU only. A CORRECTNESS REFERENCE, NOT A PERFORMANCE PATH — it is the +// numerical oracle #1189 milestone M5's CUTLASS kernel is measured against, and +// it makes no speed claim. M5 owns the CUDA arm. +void MatmulFp8BlockScaled(Queue& q, Tensor& out, const Tensor& a_fp8, const Tensor& a_scale, + const Tensor& b_fp8, const Tensor& b_scale, int block_n, + int block_k); + // RmsNormQuantFp8 (fused fp8 RMSNorm -> static per-tensor activation quant). One // HBM pass mirrors vLLM's Inductor `fused_add_rms_norm_static_fp8_quant` // (vllm/compilation/passes/fusion/rms_quant_fusion.py:124) — the RMSNorm producer diff --git a/src/vt/cpu/cpu_ops.cpp b/src/vt/cpu/cpu_ops.cpp index fcb0afd21..4768b52eb 100644 --- a/src/vt/cpu/cpu_ops.cpp +++ b/src/vt/cpu/cpu_ops.cpp @@ -628,6 +628,79 @@ void QuantFp8GroupKernel(Queue&, Tensor& out_fp8, Tensor& out_scale, const Tenso }); } +// --- Block-wise FP8 (VT-MATMUL-FP8-BLOCK-REF, #1189 M2). MatmulFp8BlockScaled +// CPU kernel: the 128x128 block-scaled fp8 GEMM, mirroring +// native_w8a8_block_matmul (tests/kernels/quant_utils.py:91-154). +// +// THE SCALES APPLY IN THE MAINLOOP, ONCE PER K-BLOCK. `part` is a SEPARATE f32 +// accumulator that is scaled and only then folded into `acc`. That separation is +// the whole point of the op: MatmulFp8CutlassKernel fifty lines below folds one +// scalar alpha AFTER the whole K reduction, which has exactly one degree of +// freedom per output element while this scheme has cdiv(K, block_k) of them. An +// epilogue-only application cannot express a per-K-block scale AT ALL. Collapsing +// `part` into `acc` here IS that epilogue form, and +// tests/vt/test_ops_matmul_fp8_block_cpu.cpp G4 is constructed so it cannot pass. +// +// The scale PRODUCT is formed first, `part * (a_s * b_s)`, mirroring +// `s = As_tiles[i] * Bs[j][i]` then `c += matmul(a, b.t()) * s` +// (quant_utils.py:150-151). The kernel that executes upstream is CUTLASS, not +// Triton (vllm/model_executor/kernels/linear/__init__.py:355-377, +// vllm/utils/deep_gemm.py:27-46), and it associates left to right -- +// `accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i)`, cutlass +// 4.5.0 sm120_mma_tma_blockwise_scaling.hpp:714-717. The two differ by at most one +// f32 ULP per K-block, and upstream's own gate is what admits it: it compares the +// CUTLASS arm against THIS reference at rel_diff < 0.001 +// (test_block_fp8.py:194-200). We mirror the reference, because this op is the +// reference port. +// +// `col / block_n` indexes the b_scale ROW by OUTPUT COLUMN, mirroring +// `offs_bsn = offs_bn // group_n` (fp8_utils.py:823) and the reference's tiling +// (quant_utils.py:131-143). For a round N that agrees with a tile counter; for a +// ragged N it does not, which is why N=576 is in the ported grid. +// +// CEIL tiling and a short final K-tile: `k1 = min(k0 + block_k, k)` +// (quant_utils.py:131-141). Upstream's wrapper asserts the ceil shapes +// (fp8_utils.py:935-936), so a ragged block is legal rather than tolerated. +// +// A CORRECTNESS REFERENCE, NOT A PERFORMANCE PATH, in the same sense as +// MatmulFp8CutlassKernel below: a naive nest that makes the block-fp8 seam +// resolvable on a CPU queue so #1189 milestones M3 and M4 can be gated without a +// GPU. It makes no speed claim. M5 owns the CUDA arm and will be measured +// against this one. Parallel over ROWS; each output row is independent and the +// reduction order inside a row is fixed, so the result does not depend on the +// thread count. +void MatmulFp8BlockScaledKernel(Queue&, Tensor& out, const Tensor& a_fp8, const Tensor& a_scale, + const Tensor& b_fp8, const Tensor& b_scale, int block_n, + int block_k) { + const int64_t m = a_fp8.shape[0], k = a_fp8.shape[1], n = b_fp8.shape[0]; + const int64_t k_tiles = (k + block_k - 1) / block_k; + const auto* ap = a_fp8.Ptr(); + const auto* bp = b_fp8.Ptr(); + const auto* asp = a_scale.Ptr(); + const auto* bsp = b_scale.Ptr(); + ForRows(m, [&](int64_t r0, int64_t r1) { + std::vector arow(static_cast(k)); + for (int64_t i = r0; i < r1; ++i) { + // Decode the A row once and reuse it across N, as MatmulNvfp4Fp4Kernel does. + for (int64_t kk = 0; kk < k; ++kk) + arow[static_cast(kk)] = Fp8ToF32(ap[i * k + kk]); + for (int64_t col = 0; col < n; ++col) { + const int64_t nb = col / block_n; // fp8_utils.py:823 + float acc = 0.0F; + for (int64_t kt = 0; kt < k_tiles; ++kt) { + const int64_t k0 = kt * block_k; + const int64_t k1 = std::min(k0 + static_cast(block_k), k); + float part = 0.0F; // the MAINLOOP register, kept separate + for (int64_t kk = k0; kk < k1; ++kk) + part += arow[static_cast(kk)] * Fp8ToF32(bp[col * k + kk]); + acc += part * (asp[i * k_tiles + kt] * bsp[nb * k_tiles + kt]); + } + StoreF32(out, i * n + col, acc); + } + } + }); +} + // MatmulFp8Cutlass CPU kernel: out[m,n] = alpha * Sum_k f8val(a[m,k])*f8val(b[n,k]), // f32 accumulate, ONE folded alpha (= input_scale*weight_scale — our recorded // deviation from upstream's two epilogue scalars, see include/vt/ops.h). @@ -3256,6 +3329,9 @@ struct Registrar { reinterpret_cast(static_cast(&QuantFp8GroupKernel))); RegisterOp(OpId::kMatmulFp8Cutlass, DeviceType::kCPU, reinterpret_cast(static_cast(&MatmulFp8CutlassKernel))); + RegisterOp(OpId::kMatmulFp8BlockScaled, DeviceType::kCPU, + reinterpret_cast( + static_cast(&MatmulFp8BlockScaledKernel))); RegisterOp(OpId::kSiluAndMul, DeviceType::kCPU, reinterpret_cast(static_cast(&SiluAndMulKernel))); RegisterOp(OpId::kGeluAndMul, DeviceType::kCPU, diff --git a/src/vt/op_provider.cpp b/src/vt/op_provider.cpp index bd91bd174..24f6432d7 100644 --- a/src/vt/op_provider.cpp +++ b/src/vt/op_provider.cpp @@ -491,6 +491,8 @@ const char* OpNameImpl(OpId op) { return "AttentionRelPos"; case OpId::kQuantFp8Group: return "QuantFp8Group"; + case OpId::kMatmulFp8BlockScaled: + return "MatmulFp8BlockScaled"; case OpId::kCount: break; } diff --git a/src/vt/ops.cpp b/src/vt/ops.cpp index 8a8f6aaa2..26b8d78d1 100644 --- a/src/vt/ops.cpp +++ b/src/vt/ops.cpp @@ -673,6 +673,50 @@ void QuantFp8Group(Queue& q, Tensor& out_fp8, Tensor& out_scale, const Tensor& x reinterpret_cast(GetOp(OpId::kQuantFp8Group, q.device.type))( q, out_fp8, out_scale, x, group_size); } +void MatmulFp8BlockScaled(Queue& q, Tensor& out, const Tensor& a_fp8, const Tensor& a_scale, + const Tensor& b_fp8, const Tensor& b_scale, int block_n, + int block_k) { + VT_CHECK(out.rank == 2 && a_fp8.rank == 2 && a_scale.rank == 2 && b_fp8.rank == 2 && + b_scale.rank == 2, + "matmul_fp8_block_scaled: out/a_fp8/a_scale/b_fp8/b_scale must be rank-2"); + // Validated BEFORE either one divides anything: `x / 0` and `x % 0` are + // undefined behaviour, so a zero must refuse rather than trap. + VT_CHECK(block_n > 0 && block_k > 0, + "matmul_fp8_block_scaled: block_n and block_k must be positive"); + const int64_t m = a_fp8.shape[0], k = a_fp8.shape[1]; + const int64_t n = b_fp8.shape[0]; + // Upstream: `assert A.shape[-1] == B.shape[-1]` (quant_utils.py:111). + VT_CHECK(b_fp8.shape[1] == k, + "matmul_fp8_block_scaled: a_fp8 [M,K] and b_fp8 [N,K] must share K"); + VT_CHECK(out.shape[0] == m && out.shape[1] == n, "matmul_fp8_block_scaled: out must be [M,N]"); + VT_CHECK(a_fp8.dtype == DType::kI8 && b_fp8.dtype == DType::kI8, + "matmul_fp8_block_scaled: a_fp8/b_fp8 must be i8 (raw fp8-e4m3fn bytes)"); + // f32, not the model dtype: upstream refuses any other scale dtype on this + // path (csrc/.../w8a8/cutlass/c3x/scaled_mm_helper.hpp:15-18) and the + // accumulator these multiply into is f32. + VT_CHECK(a_scale.dtype == DType::kF32 && b_scale.dtype == DType::kF32, + "matmul_fp8_block_scaled: a_scale/b_scale must be f32"); + VT_CHECK(out.dtype == DType::kF32 || out.dtype == DType::kBF16, + "matmul_fp8_block_scaled: out must be f32 or bf16"); + // CEIL on every tiling, so a ragged final block is legal: upstream asserts + // `triton.cdiv(N, block_n) == Bs.shape[0]` and + // `triton.cdiv(K, block_k) == Bs.shape[1]` (fp8_utils.py:935-936), and + // `triton.cdiv(A.shape[-1], block_k) == As.shape[-1]` (fp8_utils.py:930). + const int64_t k_tiles = (k + block_k - 1) / block_k; + const int64_t n_tiles = (n + block_n - 1) / block_n; + VT_CHECK(a_scale.shape[0] == m && a_scale.shape[1] == k_tiles, + "matmul_fp8_block_scaled: a_scale must be [M, cdiv(K, block_k)]"); + VT_CHECK(b_scale.shape[0] == n_tiles && b_scale.shape[1] == k_tiles, + "matmul_fp8_block_scaled: b_scale must be [cdiv(N, block_n), cdiv(K, block_k)]"); + VT_CHECK(out.IsContiguous() && a_fp8.IsContiguous() && a_scale.IsContiguous() && + b_fp8.IsContiguous() && b_scale.IsContiguous(), + "matmul_fp8_block_scaled: contiguous tensors required"); + VT_CHECK(out.device == q.device && a_fp8.device == q.device && a_scale.device == q.device && + b_fp8.device == q.device && b_scale.device == q.device, + "matmul_fp8_block_scaled: device mismatch (out/a_fp8/a_scale/b_fp8/b_scale/queue)"); + reinterpret_cast(GetOp(OpId::kMatmulFp8BlockScaled, q.device.type))( + q, out, a_fp8, a_scale, b_fp8, b_scale, block_n, block_k); +} void RmsNormQuantFp8(Queue& q, Tensor& out_fp8, Tensor* out_bf16, const Tensor& x, const Tensor& weight, const RmsNormArgs& args, Tensor* residual, float input_scale) { diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index bad02c87d..2d6367a91 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -1717,6 +1717,11 @@ vllm_cpp_add_test(test_ops_fp8_cpu vt/test_ops_fp8_cpu.cpp) # quant. CPU-gateable by construction; its CPU-vs-CUDA byte-identity arm (G6) is # CUDA-gated and reports PENDING rather than skipping where no device exists. vllm_cpp_add_test(test_ops_quant_fp8_group_cpu vt/test_ops_quant_fp8_group_cpu.cpp) +# VT-MATMUL-FP8-BLOCK-REF (#1189 M2): the 128x128 block-scaled fp8 GEMM, CPU +# reference arm — the numerical oracle every later block-fp8 kernel is measured +# against. Its G4 holds the constraint the row exists for: the scales apply in +# the MAINLOOP, once per K-block, which an epilogue-folded alpha cannot express. +vllm_cpp_add_test(test_ops_matmul_fp8_block_cpu vt/test_ops_matmul_fp8_block_cpu.cpp) # Opt-in arm: run the fp8 plan-cache byte-exact case with the cache ENABLED # (VT_FP8_PLAN_CACHE=1 -> first MatmulFp8CublasLt call builds the plan fresh, # later calls hit the cache). Proves the cached-plan GEMM is BYTE-identical to the diff --git a/tests/vt/test_ops_matmul_fp8_block_cpu.cpp b/tests/vt/test_ops_matmul_fp8_block_cpu.cpp new file mode 100644 index 000000000..480b6ecf7 --- /dev/null +++ b/tests/vt/test_ops_matmul_fp8_block_cpu.cpp @@ -0,0 +1,684 @@ +// vllm.cpp — the 128x128 block-scaled FP8 GEMM, CPU reference arm. +// +// VT-MATMUL-FP8-BLOCK-REF (.agents/specs/vt-matmul-fp8-block-ref.md), issue +// #1189 milestone M2. Pinned oracle: vLLM 5559679229bc961848b121ccdeaa8fa5d79bec98, +// asserted as the local checkout's HEAD before every anchor below was read. +// +// THE CONSTRAINT THIS FILE EXISTS TO HOLD. The scales apply in the GEMM +// MAINLOOP, once per K-block, into an f32 accumulator — not in the epilogue: +// +// accumulator = zeros(f32) +// for k_block: +// accumulator += dot(a_tile, b_tile) * a_s[m,k_block] * b_s[n_block,k_block] +// +// (vllm/model_executor/layers/quantization/utils/fp8_utils.py:826-836, the +// Triton arm spelling out the structure). Our per-tensor FP8 path folds ONE +// scalar alpha into the epilogue (`MatmulFp8CutlassKernel`, `alpha * acc` after +// the whole K reduction). An epilogue has exactly one degree of freedom per +// output element and the block scheme has cdiv(K, block_k) of them, so an +// epilogue-only application CANNOT EXPRESS a per-K-block scale at all. That is a +// correctness constraint, not an optimisation choice. G4 is the instrument for +// it and no single-alpha implementation can pass G4. +// +// WHICH IMPLEMENTATION ACTUALLY RUNS, asked because M1 found upstream's Triton +// source is not what executes and that the two arms DISAGREE in polarity. Here +// the chain lands differently — the executing kernel is CUTLASS and it AGREES +// with the reference: +// vllm/model_executor/kernels/linear/__init__.py:355-377 CUDA priority list +// vllm/utils/deep_gemm.py:27-46 DeepGEMM off for +// qwen3_5_text on family 120 +// .../scaled_mm/cutlass.py:312-326 ops.cutlass_scaled_mm(A, B.T, As, Bs.T) +// .../cutlass/scaled_mm_entry.cu:220-226 every sm>=120 -> sm120 +// .../cutlass/scaled_mm_c3x_sm120.cu:13-20 -> blockwise_sm120_fp8 +// .../c3x/scaled_mm_helper.hpp:15-18,39-55 blockwise branch: both +// scales f32 and 2-D, +// ceil_div shapes, no bias +// .../c3x/scaled_mm_blockwise_sm120_fp8_dispatch.cuh:56-58,218-235 +// ElementAccumulator=float, +// both scale pointers are +// MAINLOOP arguments +// cutlass 4.5.0 include/cutlass/gemm/collective/ +// sm120_mma_tma_blockwise_scaling.hpp:714-717 +// accum(i) += tmp_accum(i) * tCrScaleAViewAsC(i) * tCrScaleBViewAsC(i) +// which is native_w8a8_block_matmul's `c += matmul(a, b.t()) * s` +// (tests/kernels/quant_utils.py:145-151) with the two scale multiplies +// associated left to right instead of as a product. The two differ by at most +// one f32 ULP per K-block, and upstream's own gate is what admits it: it +// compares the CUTLASS arm against this very reference at rel_diff < 0.001 +// (tests/kernels/quantization/test_block_fp8.py:194-200). We mirror the +// reference's association, because this op IS the reference port. +// +// G1 the registration itself, and OpName totality. +// G2 the ported upstream case (test_block_fp8.py:123-153) over an adapted +// grid, against an INDEPENDENTLY WRITTEN double reference, carrying +// upstream's rel_diff < 0.001 verbatim AND a tighter f32 forward-error +// bound. Vacuity guard. +// G3 the ragged edges, called out separately because they are why the grid +// is what it is: N=576 (4*128+64), K=3884 (30*128+44), and both at once. +// G4 THE MAINLOOP CONSTRAINT, by exact hand-computed values plus a scale +// SWAP that no epilogue-folded alpha can distinguish. +// G5 the refusals, each by name. +// G6 the M1 seam: vt::QuantFp8Group -> vt::MatmulFp8BlockScaled end to end +// on a CPU queue. A COMPOSITION test, not a reachability claim: nothing +// in production dispatches either op yet and the spec's `## Owed` says so. +// +// THE REFERENCE IS INDEPENDENT, or the gate is a tautology. It accumulates in +// `double` with a DIFFERENT loop nest — upstream's, k-tile outer and n-block +// inner over whole tiles — while the kernel walks output elements with the +// k-tiles inside. Its fp8 decode is derived from the e4m3fn field layout +// (1 sign, 4 exponent bias 7, 3 mantissa, no infinities), not taken from any +// codec in src/. Reduction-order ULPs are exactly what the f32 forward-error +// bound admits; nothing here claims bit-exactness against a double. +// +// THE UPSTREAM GRID IS ADAPTED, and the adaptation is stated rather than +// hidden. Upstream runs itertools.product(M, N, K) with M=[1,7,8,83,4096], +// N=[128,512,576,7168,13824], K=[256,3884,4096,13824,16384] +// (test_block_fp8.py:48-55) — 125 combinations whose largest is 4096x13824x16384, +// about 9.3e11 multiply-accumulates. That is a GPU grid; a naive CPU reference +// nest, run twice, cannot execute it. Every axis VALUE below appears at least +// once and the ragged values appear separately and together; only the PAIRING +// is dropped. Parameters, dtypes, tolerance and failure criterion are preserved +// exactly. +#include + +#include +#include +#include +#include +#include +#include +#include +#include + +#include "vt/backend.h" +#include "vt/dtype.h" +#include "vt/ops.h" + +namespace { + +using vt::Device; +using vt::DeviceType; +using vt::DType; +using vt::Queue; +using vt::Tensor; + +Device Cpu() { return Device{DeviceType::kCPU, 0}; } + +Tensor MakeTensor(void* data, DType dt, Device dev, const std::vector& shape) { + Tensor t; + t.data = data; + t.dtype = dt; + t.device = dev; + t.rank = static_cast(shape.size()); + int64_t stride = 1; + for (int i = t.rank - 1; i >= 0; --i) { + t.shape[i] = shape[static_cast(i)]; + t.stride[i] = stride; + stride *= shape[static_cast(i)]; + } + return t; +} + +int64_t CDiv(int64_t a, int64_t b) { return (a + b - 1) / b; } + +// --- the independent e4m3fn decode ---------------------------------------- +// From the format: 1 sign bit, 4 exponent bits (bias 7), 3 mantissa bits, no +// infinities, 0x7F/0xFF the only NaN encodings ("fn"). Not taken from src/. +constexpr float kFp8Max = 448.0f; // quant_utils.py:27-35 finfo(e4m3fn).max +constexpr float kFp8Min = -448.0f; + +double E4m3Exact(unsigned exp_field, unsigned mant) { + if (exp_field == 0) return std::ldexp(static_cast(mant), -9); // mant / 512 + return std::ldexp(1.0 + static_cast(mant) / 8.0, static_cast(exp_field) - 7); +} + +const std::vector& DequantTable() { + static const std::vector table = [] { + std::vector v(256); + for (unsigned b = 0; b < 256; ++b) { + const double m = E4m3Exact((b >> 3) & 0xFu, b & 0x7u); + v[b] = (b & 0x80u) != 0 ? -m : m; + } + return v; + }(); + return table; +} + +// The 127 finite magnitudes 0x00..0x7E are MONOTONIC in the byte value, so the +// nearest-value encode is a binary search over their midpoints, with an exact +// midpoint resolving to the even mantissa (table index == byte value, so +// `lo & 1` IS the mantissa's low bit). +const std::vector& Magnitudes() { + static const std::vector table = [] { + std::vector v; + v.reserve(127); + for (unsigned e = 0; e <= 15; ++e) + for (unsigned m = 0; m <= 7; ++m) { + if (e == 15 && m == 7) continue; // the NaN encoding is not a value + v.push_back(E4m3Exact(e, m)); + } + return v; + }(); + return table; +} + +uint8_t EncodeRne(float r) { + const auto sign = static_cast(std::signbit(r) ? 0x80u : 0x00u); + const double a = std::fabs(static_cast(r)); + const std::vector& mag = Magnitudes(); + const auto it = std::lower_bound(mag.begin(), mag.end(), a); + size_t hi = static_cast(it - mag.begin()); + if (hi == 0) return sign; // a == 0 + if (hi >= mag.size()) return static_cast(sign | 0x7Eu); // a >= 448 + const size_t lo = hi - 1; + const double dlo = a - mag[lo], dhi = mag[hi] - a; + size_t pick = 0; + if (dlo < dhi) { + pick = lo; + } else if (dhi < dlo) { + pick = hi; + } else { + pick = (lo & 1u) == 0u ? lo : hi; + } + return static_cast(sign | static_cast(pick)); +} + +// --- operands --------------------------------------------------------------- +// Upstream builds each operand as +// A_fp32 = (rand(M,K) - 0.5) * 2 * fp8_max +// A_fp8 = A_fp32.clamp(min=fp8_min, max=fp8_max).to(fp8) +// (test_block_fp8.py:130-136), i.e. the fp8 bytes induced by a value that is +// UNIFORM on [-448, 448]. The transform is reproduced; the RNG is not torch's, +// which is a harness fact rather than a parameter — the values are random on +// both sides and nothing here depends on a particular draw. +// +// The encode is hoisted into a 2^20-entry table of uniformly spaced values +// because the grid below needs ~4e7 operand bytes and a binary search per byte, +// at the -O0 the CPU lane builds with, dominated the case's run time. The +// resulting byte distribution is upstream's to a value resolution of 8.5e-4, +// and a uniform draw lands below that magnitude with probability ~2e-6. +const std::vector& UniformByteTable() { + static const std::vector table = [] { + constexpr int kN = 1 << 20; + std::vector v(static_cast(kN)); + for (int i = 0; i < kN; ++i) { + const float u = (static_cast(i) + 0.5f) / static_cast(kN); // (0,1) + const float x = (u - 0.5f) * 2.0f * kFp8Max; + v[static_cast(i)] = EncodeRne(std::fmin(std::fmax(x, kFp8Min), kFp8Max)); + } + return v; + }(); + return table; +} + +std::vector RandomFp8(int64_t n, uint32_t seed) { + const std::vector& tab = UniformByteTable(); + std::mt19937 rng(seed); + std::uniform_int_distribution idx(0, static_cast(tab.size() - 1)); + std::vector out(static_cast(n)); + for (int64_t i = 0; i < n; ++i) out[static_cast(i)] = tab[idx(rng)]; + return out; +} + +// As = rand(M, k_tiles) * 1e-2 ; Bs = rand(n_tiles, k_tiles) * 1e-2 +// (test_block_fp8.py:143-144, factor_for_scale = 1e-2 at :127). +std::vector RandomScales(int64_t n, uint32_t seed) { + std::mt19937 rng(seed); + std::uniform_real_distribution u(0.0f, 1.0f); + std::vector out(static_cast(n)); + for (int64_t i = 0; i < n; ++i) out[static_cast(i)] = u(rng) * 1e-2f; + return out; +} + +// --- the independent reference --------------------------------------------- +// native_w8a8_block_matmul (tests/kernels/quant_utils.py:91-154), transcribed: +// k-tile outer, n-block inner, whole tiles, f32 semantics widened to double so +// that the kernel's own f32 rounding is what the bound measures. `abs_sum` +// carries the same nest and is the forward-error bound's magnitude term. +struct RefResult { + std::vector c; // [M,N] + std::vector abs_sum; // [M,N] sum over k-tiles of |partial| * |s| +}; + +RefResult RefBlockMatmul(const std::vector& ad, const std::vector& bd, + const std::vector& as, const std::vector& bs, int64_t m, + int64_t n, int64_t k, int64_t block_n, int64_t block_k) { + const int64_t n_tiles = CDiv(n, block_n); // quant_utils.py:123-124, CEIL + const int64_t k_tiles = CDiv(k, block_k); + RefResult r; + r.c.assign(static_cast(m * n), 0.0); + r.abs_sum.assign(static_cast(m * n), 0.0); + for (int64_t i = 0; i < k_tiles; ++i) { // quant_utils.py:145 + const int64_t k0 = i * block_k; + const int64_t k1 = std::min((i + 1) * block_k, k); + for (int64_t j = 0; j < n_tiles; ++j) { // quant_utils.py:146 + const int64_t n0 = j * block_n; + const int64_t n1 = std::min((j + 1) * block_n, n); + const double b_s = static_cast(bs[static_cast(j * k_tiles + i)]); + for (int64_t row = 0; row < m; ++row) { + const double a_s = static_cast(as[static_cast(row * k_tiles + i)]); + const double s = a_s * b_s; // quant_utils.py:150, PRODUCT first + for (int64_t col = n0; col < n1; ++col) { + double part = 0.0, part_abs = 0.0; + for (int64_t kk = k0; kk < k1; ++kk) { + const double p = ad[static_cast(row * k + kk)] * + bd[static_cast(col * k + kk)]; + part += p; + part_abs += std::fabs(p); + } + r.c[static_cast(row * n + col)] += part * s; // quant_utils.py:151 + r.abs_sum[static_cast(row * n + col)] += part_abs * std::fabs(s); + } + } + } + } + return r; +} + +std::vector Decode(const std::vector& bytes) { + const std::vector& tab = DequantTable(); + std::vector out(bytes.size()); + for (size_t i = 0; i < bytes.size(); ++i) out[i] = tab[bytes[i]]; + return out; +} + +double ReadOut(const std::vector& buf, DType dt, int64_t i) { + if (dt == DType::kF32) { + float v = 0.0f; + std::memcpy(&v, buf.data() + static_cast(i) * 4, sizeof(v)); + return static_cast(v); + } + uint16_t v = 0; + std::memcpy(&v, buf.data() + static_cast(i) * 2, sizeof(v)); + return static_cast(vt::BF16ToF32(v)); +} + +// One shape of the ported case. `out_dtype` is the store width; upstream runs +// bfloat16 (test_block_fp8.py:54) and the op admits f32 as well. +void RunPorted(int64_t m, int64_t n, int64_t k, DType out_dtype, uint32_t seed) { + constexpr int64_t kBlockN = 128, kBlockK = 128; // test_block_fp8.py:53 + const int64_t n_tiles = CDiv(n, kBlockN), k_tiles = CDiv(k, kBlockK); + CAPTURE(m); + CAPTURE(n); + CAPTURE(k); + + const std::vector a = RandomFp8(m * k, seed); + const std::vector b = RandomFp8(n * k, seed + 1u); + const std::vector as = RandomScales(m * k_tiles, seed + 2u); + const std::vector bs = RandomScales(n_tiles * k_tiles, seed + 3u); + + const int64_t out_elems = m * n; + std::vector out(static_cast(out_elems) * + (out_dtype == DType::kF32 ? 4u : 2u)); + + Queue q; + q.device = Cpu(); + Tensor ta = MakeTensor(const_cast(a.data()), DType::kI8, Cpu(), {m, k}); + Tensor tas = MakeTensor(const_cast(as.data()), DType::kF32, Cpu(), {m, k_tiles}); + Tensor tb = MakeTensor(const_cast(b.data()), DType::kI8, Cpu(), {n, k}); + Tensor tbs = + MakeTensor(const_cast(bs.data()), DType::kF32, Cpu(), {n_tiles, k_tiles}); + Tensor tout = MakeTensor(out.data(), out_dtype, Cpu(), {m, n}); + vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, static_cast(kBlockN), + static_cast(kBlockK)); + + const RefResult ref = RefBlockMatmul(Decode(a), Decode(b), as, bs, m, n, k, kBlockN, kBlockK); + + // (1) upstream's own criterion, verbatim (test_block_fp8.py:150-153): + // mean|out - ref| / mean|ref| < 0.001, with BOTH sides at the store width + // so that a bf16 store is common-mode exactly as it is upstream. + // (2) a tighter per-element f32 forward-error bound, because our arm IS the + // reference rather than a kernel measured against it: the f32 recursive-sum + // bound over K terms, x4 margin, plus the store's own half-ulp for bf16. + double sum_abs_diff = 0.0, sum_abs_ref = 0.0, worst_ratio = 0.0; + int64_t bad = 0, first_bad = -1, nonzero = 0; + const double eps_f32 = static_cast(std::numeric_limits::epsilon()); + for (int64_t i = 0; i < out_elems; ++i) { + const double got = ReadOut(out, out_dtype, i); + const double want = ref.c[static_cast(i)]; + const double want_stored = + out_dtype == DType::kF32 + ? static_cast(static_cast(want)) + : static_cast(vt::BF16ToF32(vt::F32ToBF16(static_cast(want)))); + sum_abs_diff += std::fabs(got - want_stored); + sum_abs_ref += std::fabs(want_stored); + const double tol = + 4.0 * static_cast(k) * eps_f32 * ref.abs_sum[static_cast(i)] + + (out_dtype == DType::kBF16 ? std::fabs(want) * 0.004 : 0.0); + const double diff = std::fabs(got - want); + if (tol > 0.0 && diff / tol > worst_ratio) worst_ratio = diff / tol; + if (!(diff <= tol)) { + if (bad == 0) { + first_bad = i; + CAPTURE(got); + CAPTURE(want); + CAPTURE(tol); + } + ++bad; + } + if (std::fabs(want) > 0.0) ++nonzero; + } + if (bad != 0) { + CAPTURE(bad); + CAPTURE(first_bad); + } + CAPTURE(worst_ratio); + CHECK(bad == 0); + const double rel_diff = sum_abs_diff / sum_abs_ref; + CAPTURE(rel_diff); + CHECK(rel_diff < 0.001); // test_block_fp8.py:153 + // VACUITY GUARD: an all-zero reference would make any implementation pass. + CHECK(nonzero == out_elems); +} + +} // namespace + +// =========================================================================== +// G1 — THE REGISTRATION ITSELF. +// +// A refusal test cannot stand in for this and that is measured, not argued: M1 +// deleted its CPU registration and every refusal case still passed, because the +// refusals live in the wrapper's validation and fire before dispatch. Only an +// OpRegistered assertion caught it (.agents/specs/vt-quant-fp8-group.md, +// "A refusal test cannot stand in for a registration test"). +TEST_CASE("G1: MatmulFp8BlockScaled is registered on the CPU backend and named") { + REQUIRE(vt::OpRegistered(vt::OpId::kMatmulFp8BlockScaled, DeviceType::kCPU)); + CHECK(std::string(vt::OpName(vt::OpId::kMatmulFp8BlockScaled)) == "MatmulFp8BlockScaled"); +} + +// =========================================================================== +// G2/G3 — the ported upstream case over the adapted grid. +// +// Every axis value of upstream's M/N/K lists appears at least once. The three +// ragged entries are G3 and are marked: N=576 is 4*128 + 64, K=3884 is +// 30*128 + 44, and one case carries both at once, which is where a floor +// division and a ceil division give different answers in different places. +TEST_CASE("G2: MatmulFp8BlockScaled matches native_w8a8_block_matmul over the ported grid") { + REQUIRE(vt::OpRegistered(vt::OpId::kMatmulFp8BlockScaled, DeviceType::kCPU)); + // M N K out seed + RunPorted(4096, 128, 256, DType::kBF16, 1101); // M=4096 + RunPorted(83, 128, 4096, DType::kBF16, 1201); // M=83, K=4096 + RunPorted(8, 7168, 256, DType::kBF16, 1301); // N=7168 + RunPorted(7, 13824, 256, DType::kBF16, 1401); // M=7, N=13824 + RunPorted(1, 512, 13824, DType::kBF16, 1501); // M=1 decode, K=13824 + RunPorted(1, 128, 16384, DType::kBF16, 1601); // K=16384 + // the f32 store arm, which upstream lists but leaves commented out at :54 + RunPorted(83, 128, 4096, DType::kF32, 1701); +} + +TEST_CASE("G3: MatmulFp8BlockScaled handles a ragged final N-block and K-block") { + REQUIRE(vt::OpRegistered(vt::OpId::kMatmulFp8BlockScaled, DeviceType::kCPU)); + RunPorted(8, 576, 4096, DType::kBF16, 2101); // N=576 = 4*128 + 64, round K + RunPorted(32, 576, 3884, DType::kBF16, 2201); // BOTH ragged: K=3884 = 30*128 + 44 + // Upstream's own dedicated ragged case, "weight.shape % 128 != 0, like in DSV3 + // kv_a_proj_with_mqa" (test_block_fp8.py:156-200). + RunPorted(32, 576, 7168, DType::kBF16, 2301); + RunPorted(32, 576, 3884, DType::kF32, 2401); +} + +// =========================================================================== +// G4 — THE MAINLOOP CONSTRAINT, and it is constructed so that no epilogue-folded +// alpha can pass it. +// +// Every value here is exact in bf16 and in f32, so these are equalities and not +// tolerances. Written as bare `==` rather than doctest's Approx: Approx carries +// a `scale` term and a strict `<`, so `epsilon(0.0)` compares `|a-b| < 0` and +// fails on values that ARE equal, while the default epsilon admits ~1.19e-5. Case A has two K-blocks whose partial products differ (128 and +// 256) and whose scales differ (0.25 and 0.5). SWAPPING the two K-block scales +// leaves every per-tensor summary of the scale tensor identical — same set, +// same sum, same product, same max — and changes the correct answer from 160 to +// 128. A kernel that reduces the scales to one epilogue alpha returns the same +// number for both and fails here. +TEST_CASE("G4: the scales apply per K-BLOCK in the mainloop, not once in the epilogue") { + REQUIRE(vt::OpRegistered(vt::OpId::kMatmulFp8BlockScaled, DeviceType::kCPU)); + constexpr uint8_t kOne = 0x38; // e4m3fn 1.0 = exponent field 7, mantissa 0 + constexpr uint8_t kTwo = 0x40; // e4m3fn 2.0 = exponent field 8, mantissa 0 + REQUIRE(DequantTable()[kOne] == 1.0); + REQUIRE(DequantTable()[kTwo] == 2.0); + + Queue q; + q.device = Cpu(); + + SUBCASE("A: swapping the two K-block scales changes the answer") { + constexpr int64_t kM = 1, kN = 1, kK = 256; // 2 K-blocks of 128, 1 N-block + std::vector a(static_cast(kM * kK), kOne); + std::vector b(static_cast(kN * kK), kOne); + for (int64_t i = 128; i < kK; ++i) b[static_cast(i)] = kTwo; // block 1: b = 2 + // partial(block 0) = 128 * 1 * 1 = 128 ; partial(block 1) = 128 * 1 * 2 = 256 + std::vector as = {1.0f, 1.0f}; + std::vector bs = {0.25f, 0.5f}; + + Tensor ta = MakeTensor(a.data(), DType::kI8, Cpu(), {kM, kK}); + Tensor tb = MakeTensor(b.data(), DType::kI8, Cpu(), {kN, kK}); + Tensor tas = MakeTensor(as.data(), DType::kF32, Cpu(), {kM, 2}); + Tensor tbs = MakeTensor(bs.data(), DType::kF32, Cpu(), {1, 2}); + + float got = 0.0f; + Tensor tout = MakeTensor(&got, DType::kF32, Cpu(), {kM, kN}); + vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, 128, 128); + // 128*0.25 + 256*0.5 = 32 + 128 + CHECK(static_cast(got) == 160.0); + + std::swap(bs[0], bs[1]); + float got_swapped = 0.0f; + Tensor tout2 = MakeTensor(&got_swapped, DType::kF32, Cpu(), {kM, kN}); + vt::MatmulFp8BlockScaled(q, tout2, ta, tas, tb, tbs, 128, 128); + // 128*0.5 + 256*0.25 = 64 + 64 + CHECK(static_cast(got_swapped) == 128.0); + // The load-bearing statement: an epilogue alpha cannot tell these apart. + CHECK(got != got_swapped); + // VACUITY GUARD. + CHECK(got != 0.0f); + CHECK(got_swapped != 0.0f); + } + + SUBCASE("B: the b-scale row is indexed by OUTPUT COLUMN / block_n, ragged block included") { + // N = 129 -> 2 N-blocks, the second of them one column wide. a = b = 1.0 + // everywhere, so every K-block partial is exactly 128 and the output is + // decided entirely by which scale pair the kernel picked. + constexpr int64_t kM = 2, kN = 129, kK = 256; + std::vector a(static_cast(kM * kK), kOne); + std::vector b(static_cast(kN * kK), kOne); + std::vector as = {1.0f, 1.0f, 2.0f, 2.0f}; // [2,2] + std::vector bs = {0.25f, 0.5f, 1.0f, 0.125f}; // [2,2] + + Tensor ta = MakeTensor(a.data(), DType::kI8, Cpu(), {kM, kK}); + Tensor tb = MakeTensor(b.data(), DType::kI8, Cpu(), {kN, kK}); + Tensor tas = MakeTensor(as.data(), DType::kF32, Cpu(), {kM, 2}); + Tensor tbs = MakeTensor(bs.data(), DType::kF32, Cpu(), {2, 2}); + + std::vector out(static_cast(kM * kN), 0.0f); + Tensor tout = MakeTensor(out.data(), DType::kF32, Cpu(), {kM, kN}); + vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, 128, 128); + + int64_t nonzero = 0; + for (int64_t row = 0; row < kM; ++row) { + const double a_row = row == 0 ? 1.0 : 2.0; + for (int64_t col = 0; col < kN; ++col) { + const size_t nb = col < 128 ? 0u : 1u; // OUTPUT COLUMN / block_n + const double want = 128.0 * (a_row * static_cast(bs[nb * 2 + 0])) + + 128.0 * (a_row * static_cast(bs[nb * 2 + 1])); + const double got = static_cast(out[static_cast(row * kN + col)]); + if (got != want) { + CAPTURE(row); + CAPTURE(col); + CAPTURE(got); + CAPTURE(want); + REQUIRE(got == want); + } + if (got != 0.0) ++nonzero; + } + } + // VACUITY GUARD. + CHECK(nonzero == kM * kN); + // Named explicitly so the ragged column is not merely covered by the loop: + // column 128 is the one-wide second N-block and must read bs row 1. + CHECK(static_cast(out[static_cast(128)]) == 144.0); + CHECK(static_cast(out[static_cast(kN + 128)]) == 288.0); + } +} + +// =========================================================================== +// G5 — the refusals, each by name. +TEST_CASE("G5: MatmulFp8BlockScaled refuses a malformed call by name") { + constexpr int64_t kM = 2, kN = 256, kK = 256; + constexpr int kBn = 128, kBk = 128; + const int64_t nt = 2, kt = 2; + std::vector a(static_cast(kM * kK), 0x38); + std::vector b(static_cast(kN * kK), 0x38); + std::vector as(static_cast(kM * kt), 1.0f); + std::vector bs(static_cast(nt * kt), 1.0f); + std::vector out(static_cast(kM * kN), 0.0f); + + Queue q; + q.device = Cpu(); + Tensor ta = MakeTensor(a.data(), DType::kI8, Cpu(), {kM, kK}); + Tensor tb = MakeTensor(b.data(), DType::kI8, Cpu(), {kN, kK}); + Tensor tas = MakeTensor(as.data(), DType::kF32, Cpu(), {kM, kt}); + Tensor tbs = MakeTensor(bs.data(), DType::kF32, Cpu(), {nt, kt}); + Tensor tout = MakeTensor(out.data(), DType::kF32, Cpu(), {kM, kN}); + // The well-formed call this case perturbs. + REQUIRE_NOTHROW(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, kBn, kBk)); + + SUBCASE("a zero or negative block size, validated BEFORE it divides anything") { + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, 0, kBk), + std::runtime_error); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, kBn, 0), + std::runtime_error); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, -128, kBk), + std::runtime_error); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs, kBn, -128), + std::runtime_error); + } + SUBCASE("a b_scale sized by FLOOR instead of cdiv, on N and on K") { + // N = 200 -> cdiv = 2, floor = 1. A floor-tiled kernel would accept this. + constexpr int64_t kNr = 200; + std::vector br(static_cast(kNr * kK), 0x38); + std::vector bsf(static_cast(1 * kt), 1.0f); + std::vector outr(static_cast(kM * kNr), 0.0f); + Tensor tbr = MakeTensor(br.data(), DType::kI8, Cpu(), {kNr, kK}); + Tensor tbsf = MakeTensor(bsf.data(), DType::kF32, Cpu(), {1, kt}); + Tensor toutr = MakeTensor(outr.data(), DType::kF32, Cpu(), {kM, kNr}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, toutr, ta, tas, tbr, tbsf, kBn, kBk), + std::runtime_error); + // K = 200 -> cdiv = 2, floor = 1, on the K axis of b_scale. + constexpr int64_t kKr = 200; + std::vector ar2(static_cast(kM * kKr), 0x38); + std::vector br2(static_cast(kN * kKr), 0x38); + std::vector as2(static_cast(kM * 2), 1.0f); + std::vector bs2(static_cast(nt * 1), 1.0f); + Tensor ta2 = MakeTensor(ar2.data(), DType::kI8, Cpu(), {kM, kKr}); + Tensor tb2 = MakeTensor(br2.data(), DType::kI8, Cpu(), {kN, kKr}); + Tensor tas2 = MakeTensor(as2.data(), DType::kF32, Cpu(), {kM, 2}); + Tensor tbs2 = MakeTensor(bs2.data(), DType::kF32, Cpu(), {nt, 1}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta2, tas2, tb2, tbs2, kBn, kBk), + std::runtime_error); + } + SUBCASE("an a_scale with the wrong group count or the wrong row count") { + std::vector narrow(static_cast(kM * 1), 1.0f); + Tensor tnarrow = MakeTensor(narrow.data(), DType::kF32, Cpu(), {kM, 1}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tnarrow, tb, tbs, kBn, kBk), + std::runtime_error); + std::vector tall(static_cast((kM + 1) * kt), 1.0f); + Tensor ttall = MakeTensor(tall.data(), DType::kF32, Cpu(), {kM + 1, kt}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, ttall, tb, tbs, kBn, kBk), + std::runtime_error); + } + SUBCASE("a scale that is not f32") { + std::vector as16(static_cast(kM * kt), 0); + Tensor tas16 = MakeTensor(as16.data(), DType::kBF16, Cpu(), {kM, kt}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas16, tb, tbs, kBn, kBk), + std::runtime_error); + std::vector bs16(static_cast(nt * kt), 0); + Tensor tbs16 = MakeTensor(bs16.data(), DType::kBF16, Cpu(), {nt, kt}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tb, tbs16, kBn, kBk), + std::runtime_error); + } + SUBCASE("a packed operand that is not i8, and a mismatched inner dimension") { + std::vector a16(static_cast(kM * kK), 0); + Tensor ta16 = MakeTensor(a16.data(), DType::kBF16, Cpu(), {kM, kK}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta16, tas, tb, tbs, kBn, kBk), + std::runtime_error); + std::vector bk(static_cast(kN * 128), 0x38); + Tensor tbk = MakeTensor(bk.data(), DType::kI8, Cpu(), {kN, 128}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, ta, tas, tbk, tbs, kBn, kBk), + std::runtime_error); + } + SUBCASE("an out that is the wrong shape, the wrong rank, or an unsupported dtype") { + std::vector small(static_cast(kM * 8), 0.0f); + Tensor tsmall = MakeTensor(small.data(), DType::kF32, Cpu(), {kM, 8}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tsmall, ta, tas, tb, tbs, kBn, kBk), + std::runtime_error); + Tensor t3d = MakeTensor(out.data(), DType::kF32, Cpu(), {1, kM, kN}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, t3d, ta, tas, tb, tbs, kBn, kBk), + std::runtime_error); + std::vector i8out(static_cast(kM * kN), 0); + Tensor ti8 = MakeTensor(i8out.data(), DType::kI8, Cpu(), {kM, kN}); + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, ti8, ta, tas, tb, tbs, kBn, kBk), + std::runtime_error); + } + SUBCASE("a non-contiguous operand") { + Tensor gappy = ta; + gappy.stride[0] = kK + 8; // a row gap: the K run is no longer contiguous + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, gappy, tas, tb, tbs, kBn, kBk), + std::runtime_error); + } + SUBCASE("a device mismatch") { + Tensor foreign = ta; + foreign.device = Device{DeviceType::kCUDA, 0}; + CHECK_THROWS_AS(vt::MatmulFp8BlockScaled(q, tout, foreign, tas, tb, tbs, kBn, kBk), + std::runtime_error); + } +} + +// =========================================================================== +// G6 — the M1 seam, end to end on a CPU queue. +// +// vt::QuantFp8Group produces exactly the a_fp8/a_scale pair this GEMM consumes, +// and this is the composition a block-FP8 linear method will run. It is NOT a +// reachability claim: nothing in production dispatches either op at this merge +// commit, milestone M4 owns the wiring, and the spec's `## Owed` records it. +TEST_CASE("G6: QuantFp8Group feeds MatmulFp8BlockScaled on a CPU queue") { + REQUIRE(vt::OpRegistered(vt::OpId::kQuantFp8Group, DeviceType::kCPU)); + REQUIRE(vt::OpRegistered(vt::OpId::kMatmulFp8BlockScaled, DeviceType::kCPU)); + constexpr int64_t kM = 5, kN = 320, kK = 256, kG = 128; + const int64_t nt = CDiv(kN, 128), kt = kK / kG; + + std::mt19937 rng(9001); + std::uniform_real_distribution ux(-3.0f, 3.0f); + std::vector x(static_cast(kM * kK)); + for (auto& v : x) v = ux(rng); + + std::vector aq(static_cast(kM * kK)); + std::vector aqs(static_cast(kM * kt)); + Queue q; + q.device = Cpu(); + Tensor tx = MakeTensor(x.data(), DType::kF32, Cpu(), {kM, kK}); + Tensor taq = MakeTensor(aq.data(), DType::kI8, Cpu(), {kM, kK}); + Tensor taqs = MakeTensor(aqs.data(), DType::kF32, Cpu(), {kM, kt}); + vt::QuantFp8Group(q, taq, taqs, tx, static_cast(kG)); + + const std::vector b = RandomFp8(kN * kK, 7001); + const std::vector bs = RandomScales(nt * kt, 7002); + std::vector out(static_cast(kM * kN), 0); + Tensor tb = MakeTensor(const_cast(b.data()), DType::kI8, Cpu(), {kN, kK}); + Tensor tbs = MakeTensor(const_cast(bs.data()), DType::kF32, Cpu(), {nt, kt}); + Tensor tout = MakeTensor(out.data(), DType::kBF16, Cpu(), {kM, kN}); + vt::MatmulFp8BlockScaled(q, tout, taq, taqs, tb, tbs, 128, static_cast(kG)); + + const RefResult ref = RefBlockMatmul(Decode(aq), Decode(b), aqs, bs, kM, kN, kK, 128, kG); + const double eps_f32 = static_cast(std::numeric_limits::epsilon()); + int64_t bad = 0, nonzero = 0; + for (int64_t i = 0; i < kM * kN; ++i) { + const double got = static_cast(vt::BF16ToF32(out[static_cast(i)])); + const double want = ref.c[static_cast(i)]; + const double tol = + 4.0 * static_cast(kK) * eps_f32 * ref.abs_sum[static_cast(i)] + + std::fabs(want) * 0.004; + if (!(std::fabs(got - want) <= tol)) ++bad; + if (std::fabs(want) > 0.0) ++nonzero; + } + CHECK(bad == 0); + // VACUITY GUARD. + CHECK(nonzero == kM * kN); +} From d6c3a092caed03b5878b649f206e291b188172c8 Mon Sep 17 00:00:00 2001 From: Ettore Di Giacinto Date: Tue, 18 Aug 2026 11:02:26 +0000 Subject: [PATCH 3/3] record(VT-MATMUL-FP8-BLOCK-REF): the ragged shapes are the gate -- two floor-vs-ceil mutations left the round-shape grid entirely green Fills in the spec's `## Evidence` with what the row measured rather than what it intended, taken on the merged tree with `origin/main` at `2d26da5a1`. RED: with the test present and the implementation reversed out of the four source files -- `git diff --stat` printing `4 files changed, 200 deletions(-)`, so the revert is proven to have landed rather than assumed -- the focused build fails `compile_rc=1` with 29 errors, 23 `'MatmulFp8BlockScaled' is not a member of 'vt'` and 6 `'kMatmulFp8BlockScaled' is not a member of 'vt::OpId'`. GREEN: 6 cases, 80 assertions, 0 failed. The per-block table sums to 80, so no block is silently empty, and each block was filtered through a `-tc` prefix containing no comma because doctest splits that flag on commas and a name carrying one reports `0 cases ran` under a `SUCCESS!` banner. The finding worth keeping is about the grid rather than the kernel. Two separate floor-versus-ceil mutations -- `k_tiles = k / block_k` in the kernel, and a b_scale row index of `min(col / block_n, n / block_n - 1)` -- both compiled, both ran, and both left **G2 entirely green**, because every N and every K in the round-shape block is a multiple of 128. Only G3 caught either one. That is the argument for carrying upstream's non-round shapes, `N=576` (`4*128 + 64`) and `K=3884` (`30*128 + 44`), separately and together, instead of the shapes the target checkpoint happens to have. The epilogue mutation is recorded because it is the row's whole claim: folding the per-block scales into one alpha compiles, runs, and returns the identical number for a scale tensor and for that same tensor with its two K-block entries swapped. G4 is the only block that can see it, and it does. Two instrument results ride along. Deleting the block-size positivity check took SIGFPE on integer division, and doctest had already printed `assertions: 52 | 52 passed | 0 failed` before the process died, so the run's exit status is the verdict and its summary is not. And a direct run mid-session reproduced the zeroed-store mutation's exact signature against clean sources: the binary was the one the mutation harness had last linked, and the whole-tree build that followed had not reached the test target. That false red is recorded with its control -- require `ninja: no work to do` and an executable newer than every source it links -- because it cost a cycle and will cost the next reader one too. FOLLOWING_AGENTS_PROTOCOL Following-Agents-Protocol: true AI-Assisted: true Assisted-by: AGENT:claude-opus-5-1m [claude-code] --- .agents/specs/vt-matmul-fp8-block-ref.md | 100 +++++++++++++++++++++-- 1 file changed, 95 insertions(+), 5 deletions(-) diff --git a/.agents/specs/vt-matmul-fp8-block-ref.md b/.agents/specs/vt-matmul-fp8-block-ref.md index 1d5195e6c..32822f14b 100644 --- a/.agents/specs/vt-matmul-fp8-block-ref.md +++ b/.agents/specs/vt-matmul-fp8-block-ref.md @@ -289,8 +289,98 @@ scoped so that it does not. ## Evidence -Recorded below at completion: the RED before the implementation existed, the -GREEN after, the per-block gate counts, and the mutation table. Every mutation -prints `compile_rc` **and** `git diff --stat`, because a mutation that fails to -build and a mutation that never applied both read as a passing test, and M1 hit -each of those once. +Taken on the merged tree, `origin/main` at `2d26da5a1`. + +**RED.** With the test present and the implementation reversed out — the reverse +of the implementation diff applied to the four source files, `git diff --stat` +printing `4 files changed, 200 deletions(-)` so the revert is proven to have +landed — the focused build fails, `compile_rc=1`, with 29 errors: 23 +``'MatmulFp8BlockScaled' is not a member of 'vt'`` and 6 +``'kMatmulFp8BlockScaled' is not a member of 'vt::OpId'``. Restoring the four +files returns `git status --porcelain` to empty and the build to green. + +**GREEN.** `test_ops_matmul_fp8_block_cpu` reports 6 cases, 80 assertions, 0 +failed, in 3.5 s under `ctest`. Per block, each run through a `-tc` prefix +filter that contains no comma, because doctest splits `-tc` on commas and a +name that contains one yields `0 cases ran` under a `SUCCESS!` banner: + +| Block | Cases | Assertions | +|---|---:|---:| +| G1 registration and name | 1 | 2 | +| G2 the ported grid | 1 | 22 | +| G3 the ragged edges | 1 | 13 | +| G4 the mainloop constraint | 1 | 14 | +| G5 the refusals | 1 | 25 | +| G6 the M1 seam | 1 | 4 | +| **sum** | **6** | **80** | + +The buckets sum to the whole-run count, so no block is silently empty and no +filter selected nothing. + +The other declared gates on the same tree: `test_op_provider` passed, +`test_ops_quant_fp8_group_cpu` passed (M1 unchanged), +`test_ops_fp8_cpu` passed (the per-tensor sibling unchanged). +`scripts/agent-preflight.sh --fail-on-skip` reports **All gates green** with 80 +`ok` results, no `FAIL` and no `SKIP` before the verdict line. +`test_cpu_x86_llamacpp_floor` (#618) passed on this run rather than reporting +`NO_QUIET_WINDOW`, so this row needed no pristine-baseline reproduction. + +### The ragged edge, specifically + +`N=576` is `4*128 + 64` and `K=3884` is `30*128 + 44`, the two non-round shapes +from upstream's own grid (`test_block_fp8.py:49-50`). G3 runs each separately +and both together, plus upstream's dedicated DSV3 `kv_a_proj_with_mqa` case +`M=32, N=576, K=7168` (`test_block_fp8.py:156-200`). This is load-bearing and +measured, not asserted: mutation M2 (floor `k_tiles`) and mutation M3 (a floor +N-block index) each leave **G2 entirely green** and fail only G3, because every +`N` and every `K` in G2 is a multiple of 128. A grid of round shapes passes +while being wrong. + +### Mutation results + +Every mutation printed `compile_rc` **and** `git diff --stat` before the run, +because a mutation that fails to build and a mutation that never applied both +read as a passing test, and M1 of #1189 hit each of those once. Each was +restored with `git checkout --` and the restore verified. + +| Mutation | `compile_rc` | Result | +|---|---|---| +| the scales folded into ONE epilogue alpha, `part` collapsed into `acc` | 0 | 4 of 6 cases fail, 26 of 80 assertions. G4's `got != got_swapped` fires, which is the epilogue signature exactly | +| `k_tiles = k / block_k`, floor instead of cdiv, in the kernel | 0 | **only G3 fails**, 4 assertions. G2 stays green: every K in it is a multiple of 128 | +| the b_scale row index forced to 0 | **1** | proves nothing: `-Werror=unused-parameter` on `block_n` | +| the same defect with the parameter kept live: `nb = min(col / block_n, n / block_n - 1)`, a floor tile count that makes a short final N-block reuse the previous scale row | 0 | 3 of 6 cases fail, 13 assertions: G3, G4, G6. **G2 stays green**, for the same reason | +| `a_scale` dropped from the scale product | 0 | 4 of 6 cases fail, 24 assertions | +| the CPU registration deleted, the kernel symbol kept live | 0 | all 6 cases fail after 7 assertions. G1's `REQUIRE(OpRegistered(...))` is what catches it | +| `n_tiles = n / block_n`, floor, in the wrapper's validation | 0 | 4 of 6 cases fail. G5's floor-sized `b_scale` refusal stops firing and G3's well-formed ragged calls are refused instead | +| the `block_n > 0 && block_k > 0` check removed | 0 | **SIGFPE, core dumped, `run_rc=136`**. The check converts undefined behaviour into a named refusal | +| the `a_scale` shape check widened to a rank check | 0 | G5's narrow and tall `a_scale` refusals both stop firing, 2 assertions | +| the store zeroed, `acc * 0.0F` | 0 | 4 of 6 cases fail, 29 assertions, including G4's two vacuity guards — the guards are live | + +Three results are worth keeping. + +**The epilogue mutation is the row's whole claim, and it is now measured.** +Folding the per-block scales into a single alpha compiles, runs, and produces the +identical number for a scale tensor and for that same tensor with its two +K-block entries swapped. G4 is the only block that can see that, and it does. + +**A grid of round shapes is blind to the ragged bug in both directions.** Two +different floor-vs-ceil mutations — one on the K tiling in the kernel, one on the +N-block index — left G2 completely green and were caught only by G3. That is the +argument for keeping `N=576` and `K=3884` in the grid rather than the shapes the +target checkpoint happens to use, which are all multiples of 128. + +**A refusal removed can crash rather than report.** With the block-size +positivity check deleted, the process took SIGFPE on integer division and doctest +had already printed `assertions: 52 | 52 passed | 0 failed` before it died. The +run's exit status is the verdict; its printed summary is not. + +### A false red, recorded because it cost a cycle + +Mid-run, a direct invocation of the test binary reproduced the zeroed-store +mutation's exact signature — 4 failed cases, 29 failed assertions — against a +tree whose sources were clean. The binary was the one the mutation harness had +last linked; the whole-tree build that followed had not reached the test target +before the run. Nothing was wrong with the code. The control is the one +`.agents/verification.md` already names: rebuild explicitly, require +`ninja: no work to do`, and check that the executable's mtime is newer than +every source it links, before believing any verdict about it.