From 3e5bb715e8bf27508a5825afca83473bebfd90d9 Mon Sep 17 00:00:00 2001 From: Pekka Heikura Date: Tue, 16 Jun 2026 14:22:40 +0300 Subject: [PATCH 1/2] =?UTF-8?q?perf(cuda):=20Q6=5FK=20decode=20MMQ=20via?= =?UTF-8?q?=20SoA=20reader-port=20(free=20AoS)=20=E2=80=94=20+17%=20N=3D8?= =?UTF-8?q?=20(#204,=20#203)?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Routes the Q6_K trunk weights (lm-head + ffn_down-half) through the int8 tensor-core decode-MMQ tile, the last big batched-decode chunk (~8.3 ms of the ~25 ms N=8 step on the bit-exact WS-sw matvec). Measured on Qwen3-8B Q4_K_M @ 4070 Ti (aggregate decode t/s, vs the AoS/WS-sw baseline): N=6 269.0 -> 304.0 (+13.0%) N=7 305.1 -> 346.9 (+13.7%) N=8 327.2 -> 382.1 (+16.8%, reproduced 324.7 -> 379.6) peak VRAM 11394 / 12282 MiB — no spill. The decode MMQ needs coalesced staging (a SoA weight layout), which on a 12 GB card is fatal if kept as a *companion* beside the interleaved AoS weight (it ~2.3x the Q6_K weight VRAM and spilled at N=8: 325 -> 169 t/s). Earlier zero-VRAM AoS-direct staging avoided the spill but regressed -7..-11% (uncoalesced 210-B block stride). The fix is the reader-port: repack each 2-D Q6_K matmul weight to the int8 SoA layout ONCE at upload and FREE the AoS copy (like Q4_K's RepackQ4KSoa), so SoA is the only copy (~+0.4 GB net, fits N=8). To free the AoS, every Q6_K reader now consumes SoA: new _soa variants of llm_matvec_q6k / _n2 / _gemm_n / _ws / _ws_sw and llm_dequant_q6k_to_f16, each BIT-IDENTICAL to its AoS counterpart (SoA stores (q6-32) int8 + verbatim scales + d, so d*scale*(q6-32) and the reduction order are unchanged). Dispatch picks the _soa kernel when the weight handle is SoA-repacked. The token embedding table (read by embed_lookup, untied from output.weight) is excluded from the repack and stays AoS. SHARPI_Q6K_SOA=0 reverts to AoS everywhere (kill-switch); SHARPI_Q6K_DECODE_MMQ=0 keeps SoA readers but disables just the decode MMQ tile. Bit-exactness verified: CudaGemmQ6KTests.Q6KSoaReaders_AreBitIdenticalToAos + MatMulBatchedGemm_Q6K (prefill), CudaDecodeMmqTests Q6_K vs SimdKernels.DotQ6K, all 8 CudaBatchForwardMultiTests, and all 5 CudaSpecBatchVerifyTests incl. the 48-token single-user SpecDecode_GreedyParity_E2E (the ported single-user matvec stays byte-stable). Build clean (TreatWarningsAsErrors). Co-Authored-By: Claude Opus 4.8 --- src/SharpInference.Cuda/CudaBackend.cs | 134 +++++++- src/SharpInference.Cuda/CudaTextKernels.cs | 303 ++++++++++++++++++ src/SharpInference.Cuda/CudaWsKernels.cs | 274 +++++++++++++++- src/SharpInference.Engine/CudaForwardPass.cs | 28 ++ .../CudaDecodeMmqTests.cs | 91 ++++++ .../CudaGemmQ6KTests.cs | 113 +++++++ 6 files changed, 924 insertions(+), 19 deletions(-) diff --git a/src/SharpInference.Cuda/CudaBackend.cs b/src/SharpInference.Cuda/CudaBackend.cs index bce4874..7ac4265 100644 --- a/src/SharpInference.Cuda/CudaBackend.cs +++ b/src/SharpInference.Cuda/CudaBackend.cs @@ -143,8 +143,10 @@ public sealed unsafe class CudaBackend : IComputeBackend, IImageOpsBackend, IDis private nint _matvecQ4KKernel; private nint _matvecQ4KSoaKernel; // #156: scale-pre-unpacked SoA decode matvec private nint _q4kRepackSoaKernel; // #156: one-time Q4_K → SoA repack + private nint _q6kRepackSoaKernel; // #204: one-time Q6_K → SoA repack for decode MMQ private nint _matvecQ5KKernel; private nint _matvecQ6KKernel; + private nint _matvecQ6KSoaKernel; // #204: bit-identical Q6_K decode matvec over the SoA layout // Q4_0 matvec (issue #124, Gemma 4 12B QAT): keeps the q4_0 weights packed on // the GPU. Without it q4_0 falls to the F32-dequant upload (~4× VRAM — a 7 GB // model would need ~28 GB, defeating full offload). 8 rows/block × 32 thr/row. @@ -172,6 +174,7 @@ public sealed unsafe class CudaBackend : IComputeBackend, IImageOpsBackend, IDis private nint _matvecQ4KN2SoaKernel; // #156: N=2 over scale-pre-unpacked SoA weight private nint _matvecQ5KN2Kernel; private nint _matvecQ6KN2Kernel; + private nint _matvecQ6KN2SoaKernel; // #204: N=2 over the Q6_K SoA layout // Issue #111: batched GEMM-N variants — one weight matrix, N input vectors, // N output rows in a single launch. Each (row, token) runs the identical // per-row reduction as the GEMV so results are bit-identical to N sequential @@ -182,6 +185,7 @@ public sealed unsafe class CudaBackend : IComputeBackend, IImageOpsBackend, IDis private nint _matvecQ4KGemmNSoaKernel; // #156: GEMM-N over scale-pre-unpacked SoA weight private nint _matvecQ5KGemmNKernel; private nint _matvecQ6KGemmNKernel; + private nint _matvecQ6KGemmNSoaKernel; // #204: GEMM-N over the Q6_K SoA layout private nint _matvecQ80GemmNKernel; // Issue #194: weight-stationary small-N batched-decode matvecs — token loop inside // the block so each weight read is amortized across the batch. One handle per @@ -195,12 +199,25 @@ public sealed unsafe class CudaBackend : IComputeBackend, IImageOpsBackend, IDis private readonly nint[] _matvecQ80WsSoaKernels = new nint[CudaWsKernels.Variants.Length]; // #201: scale-word Q6_K WS variant + Q4_K decode MMQ (BN=16 int8 mma). private readonly nint[] _matvecQ6KWsSwKernels = new nint[CudaWsKernels.Variants.Length]; + // #204: SoA Q6_K WS variant (bit-identical to both AoS WS variants; used when the + // Q6_K weight has been repacked to SoA, which is now always for 2-D trunk weights). + private readonly nint[] _matvecQ6KWsSoaKernels = new nint[CudaWsKernels.Variants.Length]; private nint _mmqQ4kSoaActsN16Kernel; private nint _mmqQ4kSoaActsN16Bm32Kernel; // #205: BM=32 tile for grid-starved low-row shapes + // #204: Q6_K decode MMQ tiles (BN=16 int8 mma, BM=64 default + BM=32 for low-row shapes). + // RepackQ6KSoa frees the interleaved AoS weight (the SoA buffer is the only copy), so every + // Q6_K reader — including this decode-MMQ tile — reads the SoA layout in place; there is no + // AoS-direct decode-MMQ variant any more. + private nint _mmqQ6kSoaActsN16Kernel; + private nint _mmqQ6kSoaActsN16Bm32Kernel; // #205 kill-switch: SHARPI_DECODE_MMQ_BM32=0 forces the BM=64 decode-MMQ tile for all // shapes (BM=32 is default-on for grid-starved low-row shapes; output is bit-identical). private readonly bool _decodeMmqBm32Enabled = Environment.GetEnvironmentVariable("SHARPI_DECODE_MMQ_BM32") != "0"; + // #204 kill-switch: SHARPI_Q6K_DECODE_MMQ=0 disables the Q6_K decode-MMQ tile (the + // Q6_K trunk shapes fall back to the bit-exact weight-stationary matvec). Default on. + private readonly bool _q6kDecodeMmqEnabled = + Environment.GetEnvironmentVariable("SHARPI_Q6K_DECODE_MMQ") != "0"; // Issue #141: compute-bound prefill GEMM — dequant Q8_0 weight + convert // activations to fp16, then one cublasGemmEx (weight read once per batch). private nint _dequantQ80F16Kernel; @@ -211,6 +228,7 @@ public sealed unsafe class CudaBackend : IComputeBackend, IImageOpsBackend, IDis // attn_v in Q6_K; without this the Q6_K trunk matmuls fell to the per-token GEMM-N // matvec (weight re-streamed once/token), the dominant large-N prefill cost. private nint _dequantQ6KF16Kernel; + private nint _dequantQ6KF16SoaKernel; // #204: dequant over the Q6_K SoA layout private nint _dequantQ5KF16Kernel; // #162: same path for Q5_K_M mixes private nint _dequantQ40F16Kernel; // #124: Q4_0 weight → fp16 (Gemma 4 12B QAT) private nint _f32ToF16Kernel; @@ -784,6 +802,7 @@ public void Free(Tensor tensor) // the set doesn't grow across model load/free cycles. _soaHandles.TryRemove(tensor.Handle, out _); _soaQ4kHandles.TryRemove(tensor.Handle, out _); // #156 + _soaQ6kHandles.TryRemove(tensor.Handle, out _); // #204 (the repacked Q6_K SoA weight) _soaQ40Handles.TryRemove(tensor.Handle, out _); // #124/#173 if (_viewHandles.TryRemove(tensor.Handle, out _)) @@ -1873,11 +1892,12 @@ public void MatMul(Tensor output, Tensor matrix, Tensor vector, DType weightDTyp bool soa = weightDType == DType.Q8_0 && _soaHandles.ContainsKey(matrix.Handle); bool soaQ40 = weightDType == DType.Q4_0 && _soaQ40Handles.ContainsKey(matrix.Handle); + bool soaQ6k = weightDType == DType.Q6_K && _soaQ6kHandles.ContainsKey(matrix.Handle); nint kernel = weightDType switch { DType.Q4_0 => soaQ40 ? _matvecQ40SoaKernel : _matvecQ40Kernel, DType.Q5_K => _matvecQ5KKernel, - DType.Q6_K => _matvecQ6KKernel, + DType.Q6_K => soaQ6k ? _matvecQ6KSoaKernel : _matvecQ6KKernel, DType.Q8_0 => soa ? _matvecQ80SoaKernel : _matvecQ80Kernel, DType.Float32 => _matvecF32Kernel, _ => throw new NotSupportedException($"CUDA MatMul: weight dtype {weightDType} not supported (expected Q4_0, Q4_K, Q5_K, Q6_K, Q8_0, or Float32)."), @@ -1941,10 +1961,11 @@ public void MatMulN2(Tensor outputA, Tensor outputB, (nint)(&pRows), (nint)(&pCols) }; + bool soaQ6k = weightDType == DType.Q6_K && _soaQ6kHandles.ContainsKey(matrix.Handle); nint kernel = weightDType switch { DType.Q5_K => _matvecQ5KN2Kernel, - DType.Q6_K => _matvecQ6KN2Kernel, + DType.Q6_K => soaQ6k ? _matvecQ6KN2SoaKernel : _matvecQ6KN2Kernel, DType.Float32 => _matvecF32N2Kernel, _ => throw new NotSupportedException( $"CUDA MatMulN2: weight dtype {weightDType} not supported (expected Q4_K, Q5_K, Q6_K, or Float32)."), @@ -2018,9 +2039,10 @@ public void MatMulBatched(Tensor outputAll, Tensor matrix, Tensor inputAll, // All take F32 input; the Q5_K/Q6_K/Q8_0 kernels decode the weight per // element. Same (rows+7)/8 × nTok geometry across all four. bool soa = weightDType == DType.Q8_0 && _soaHandles.ContainsKey(matrix.Handle); + bool soaQ6k = weightDType == DType.Q6_K && _soaQ6kHandles.ContainsKey(matrix.Handle); nint kernel = weightDType switch { - DType.Q6_K => _matvecQ6KGemmNKernel, + DType.Q6_K => soaQ6k ? _matvecQ6KGemmNSoaKernel : _matvecQ6KGemmNKernel, DType.Q5_K => _matvecQ5KGemmNKernel, DType.Q8_0 => soa ? _matvecQ80GemmNSoaKernel : _matvecQ80GemmNKernel, _ => _matvecF32GemmNKernel, @@ -2115,12 +2137,15 @@ weightDType is not (DType.Q4_K or DType.Q5_K or DType.Q6_K or DType.Q8_0 or DTyp } // F32-input kernels: same (rows+7)/8 × 256-thread geometry as the GEMM-N - // dispatch, minus the token grid dimension. #201: Q6_K defaults to the - // scale-word variant (bit-identical, word-loaded scale/d tail). + // dispatch, minus the token grid dimension. #204: a SoA-repacked Q6_K weight + // routes to the bit-identical SoA WS reader; otherwise #201's scale-word variant + // (default) or the plain #194 AoS WS variant. bool q80Soa = weightDType == DType.Q8_0 && _soaHandles.ContainsKey(matrix.Handle); + bool q6kSoa = weightDType == DType.Q6_K && _soaQ6kHandles.ContainsKey(matrix.Handle); nint kernel = weightDType switch { - DType.Q6_K => WsV2Enabled ? _matvecQ6KWsSwKernels[variant] : _matvecQ6KWsKernels[variant], + DType.Q6_K => q6kSoa ? _matvecQ6KWsSoaKernels[variant] + : WsV2Enabled ? _matvecQ6KWsSwKernels[variant] : _matvecQ6KWsKernels[variant], DType.Q5_K => _matvecQ5KWsKernels[variant], DType.Q8_0 => q80Soa ? _matvecQ80WsSoaKernels[variant] : _matvecQ80WsKernels[variant], _ => _matvecF32WsKernels[variant], @@ -2219,13 +2244,26 @@ public void MatMulBatchedDecodeMmq(Tensor outputAll, Tensor matrix, Tensor input int rows = (int)(outputAll.ElementCount / nTok); int cols = (int)(inputAll.ElementCount / nTok); - if (weightDType != DType.Q4_K || !_soaQ4kHandles.ContainsKey(matrix.Handle) - || (cols & 0xff) != 0 || nTok < 5 || rows < 2048) + // Q4_K (#201/#205) and Q6_K (#204) both have a SoA-repacked-weight decode-MMQ tile; + // a mixed-dtype trunk (Qwen3-8B Q4_K_M carries Q6_K ffn_down-half/attn_v/lm-head) + // routes per tensor. Same eligibility floor as Q4_K: cols % 256, N ≥ 5, rows ≥ 2048 + // (below 2048 the (rows/64)-block grid starves — attn_v rows=1024 stays WS). + bool q4kEligible = weightDType == DType.Q4_K && _soaQ4kHandles.ContainsKey(matrix.Handle); + // #204: the Q6_K weight is now repacked to SoA in place (RepackQ6KSoa frees the AoS), so + // the decode-MMQ tile reads `matrix` directly — same as Q4_K. Eligibility = the weight is + // SoA-registered + the tile compiled + the kill-switch is on (SHARPI_Q6K_DECODE_MMQ). + // Same floor as Q4_K: cols % 256, N ≥ 5, rows ≥ 2048 (below 2048 the (rows/64)-block grid + // starves — attn_v rows=1024 stays WS). + bool q6kEligible = weightDType == DType.Q6_K && _q6kDecodeMmqEnabled + && _mmqQ6kSoaActsN16Kernel != nint.Zero + && _soaQ6kHandles.ContainsKey(matrix.Handle); + if ((!q4kEligible && !q6kEligible) || (cols & 0xff) != 0 || nTok < 5 || rows < 2048) { MatMulBatchedWeightStationary(outputAll, matrix, inputAll, nTok, weightDType); return; } + // Both Q4_K and Q6_K read the in-place SoA weight directly (no companion). nint wPtr = GetDevPtr(matrix); nint xPtr = GetDevPtr(inputAll); nint yPtr = GetDevPtr(outputAll); @@ -2269,14 +2307,20 @@ public void MatMulBatchedDecodeMmq(Tensor outputAll, Tensor matrix, Tensor input // 151936) keep BM=64 — they already fill the grid and BM=32 would double the // per-block activation re-staging. Output is bit-identical between the two tiles. // _smCount == 0 (attribute query failed) disables the BM=32 route (keeps BM=64). - bool useBm32 = _decodeMmqBm32Enabled && _mmqQ4kSoaActsN16Bm32Kernel != nint.Zero + // #204: Q6_K uses its SoA decode-MMQ tile (the SoA weight is the only copy); Q4_K + // uses its SoA tile. Both are bit-identical between BM=64 / BM=32. + nint bm64 = q6kEligible ? _mmqQ6kSoaActsN16Kernel : _mmqQ4kSoaActsN16Kernel; + nint bm32 = q6kEligible ? _mmqQ6kSoaActsN16Bm32Kernel : _mmqQ4kSoaActsN16Bm32Kernel; + bool useBm32 = _decodeMmqBm32Enabled && bm32 != nint.Zero && _smCount > 0 && (rows + 63) / 64 < 2 * _smCount; - nint kernel = useBm32 ? _mmqQ4kSoaActsN16Bm32Kernel : _mmqQ4kSoaActsN16Kernel; + nint kernel = useBm32 ? bm32 : bm64; uint gx = useBm32 ? (uint)((rows + 31) / 32) : (uint)((rows + 63) / 64); uint block = useBm32 ? 128u : 256u; int rm = NvrtcInterop.LaunchKernel(kernel, gx, gy, 1, block, 1, 1, 0, _stream, args, null); - if (rm != 0) throw new InvalidOperationException($"cuLaunchKernel(mmq_q4k_soa_acts_n16{(useBm32 ? "_bm32" : "")}) failed: {rm}"); + string kname = q6kEligible ? "q6k_soa" : "q4k_soa"; + if (rm != 0) throw new InvalidOperationException( + $"cuLaunchKernel(mmq_{kname}_acts_n16{(useBm32 ? "_bm32" : "")}) failed: {rm}"); } } @@ -2344,7 +2388,7 @@ public void MatMulBatchedGemm(Tensor outputAll, Tensor matrix, Tensor inputAll, nint dqKern = weightDType switch { DType.Q4_K => _soaQ4kHandles.ContainsKey(matrix.Handle) ? _dequantQ4KF16SoaKernel : _dequantQ4KF16Kernel, - DType.Q6_K => _dequantQ6KF16Kernel, // #162 + DType.Q6_K => _soaQ6kHandles.ContainsKey(matrix.Handle) ? _dequantQ6KF16SoaKernel : _dequantQ6KF16Kernel, // #162/#204 DType.Q5_K => _dequantQ5KF16Kernel, // #162 DType.Q4_0 => _soaQ40Handles.ContainsKey(matrix.Handle) ? _dequantQ40F16SoaKernel : _dequantQ40F16Kernel, // #124/#173 _ => _soaHandles.ContainsKey(matrix.Handle) ? _dequantQ80F16SoaKernel : _dequantQ80F16Kernel, @@ -2591,6 +2635,15 @@ public void MatMulBatchedMmqSoa(Tensor outputAll, Tensor soaWeight, Tensor input /// matvec auto-routes these to llm_matvec_q4k_soa. private readonly ConcurrentDictionary _soaQ4kHandles = new(); + /// Issue #204: device handles of Q6_K weights repacked into the + /// scale-pre-unpacked SoA layout ([Q (q6−32) int8][S int8 scales][D fp16 d], see + /// ). Like Q4_K (#156), FREES the + /// interleaved AoS weight — the SoA buffer is the only copy — so EVERY Q6_K reader + /// (single-token / N=2 / GEMM-N matvec, WS matvec, prefill GEMM dequant, and the + /// decode-MMQ tile) auto-routes to its *_soa + /// kernel for these handles. + private readonly ConcurrentDictionary _soaQ6kHandles = new(); + /// Issue #124/#173: device handles of Q4_0 weights repacked into the SoA /// layout (see ). The MMQ / decode dp4a / fp32 matvec / /// GEMM-fallback dequant all auto-route to the aligned-load SoA kernels. @@ -2657,6 +2710,45 @@ public Tensor RepackQ4KSoa(Tensor src, int rows, int cols) return dst; } + /// + /// Issue #204: allocate a new buffer and repack the interleaved Q6_K weight + /// [rows × nb × 210 B] into the scale-pre-unpacked SoA layout + /// [Q rows*nb*256][S rows*nb*16][D rows*nb*4] (see llm_q6k_repack_soa), FREE + /// , and mark the new handle so every Q6_K reader routes to its + /// *_soa kernel. The Q region stores the signed int8 (q6 − 32) per natural + /// element (the matvec's pre-multiply weight); S the 16 int8 scales verbatim; D the fp16 d. + /// + /// Mirrors : the interleaved weight is freed so the SoA + /// buffer is the ONLY copy (net only ~+0.4 GB over the 210 B/super-block AoS — the Q region + /// grows from 192 to 256 B and S/D round up). Every reader is SoA-aware + /// (llm_matvec_q6k_soa, ..._n2_soa, ..._gemm_n_soa, + /// ..._ws_soa, llm_dequant_q6k_to_f16_soa, and the decode-MMQ tile + /// llm_mmq_q6k_soa_acts_n16), all bit-identical to their AoS counterparts. + /// + public Tensor RepackQ6KSoa(Tensor src, int rows, int cols) + { + EnsureImageKernels(); + if (!_imageKernelsAvailable) + throw new NotSupportedException("NVRTC kernels are not available on this system."); + if ((cols & 0xff) != 0) + throw new InvalidOperationException($"RepackQ6KSoa requires cols % 256 == 0 (got {cols})."); + + long nb = cols / 256; + long totalSub = (long)rows * nb; + long bytes = totalSub * (256L + 16L + 4L); // Q + S + D regions + var dst = AllocateRawBytes(bytes, DType.Q6_K, exact: true); + nint sPtr = GetDevPtr(src), dPtr = GetDevPtr(dst); + int pRows = rows, pCols = cols; + nint* args = stackalloc nint[4] { (nint)(&sPtr), (nint)(&dPtr), (nint)(&pRows), (nint)(&pCols) }; + uint grid = (uint)((totalSub + 255) / 256); + int r = NvrtcInterop.LaunchKernel(_q6kRepackSoaKernel, grid, 1, 1, 256, 1, 1, 0, _stream, args, null); + if (r != 0) throw new InvalidOperationException($"cuLaunchKernel(q6k_repack_soa) failed: {r}"); + Synchronize(); + Free(src); + _soaQ6kHandles[dst.Handle] = 0; + return dst; + } + /// /// Issue #124/#173: allocate a new buffer the same size as the interleaved Q4_0 /// [rows×nb×18 B], repack it into the SoA layout @@ -6201,9 +6293,12 @@ private void ForceEagerJit() _embedLookupQ80Kernel, _embedLookupQ80BatchedKernel, _dequantRowsQ80Kernel, _dequantRowsQ6KKernel, _matvecF32Kernel, _matvecQ40Kernel, _matvecQ4KKernel, _matvecQ5KKernel, _matvecQ6KKernel, + _matvecQ6KSoaKernel, // #204 _matvecQ80Kernel, _matvecF32N2Kernel, _matvecQ4KN2Kernel, _matvecQ5KN2Kernel, _matvecQ6KN2Kernel, + _matvecQ6KN2SoaKernel, // #204 _matvecF32GemmNKernel, _matvecQ4KGemmNKernel, _matvecQ5KGemmNKernel, _matvecQ6KGemmNKernel, + _matvecQ6KGemmNSoaKernel, // #204 _matvecQ80GemmNKernel, _mmqQ80Kernel, _mmqQ80SoaKernel, _mmqQ4kKernel, _mmqQ4kSoaKernel, _matvecQ80Dp4aSoaKernel, _q80RepackSoaKernel, _matvecQ80SoaKernel, _matvecQ80GemmNSoaKernel, _dequantQ80F16SoaKernel, @@ -6211,7 +6306,7 @@ private void ForceEagerJit() // AoS dequant (all were missing — eager-JIT them so first decode pays no stutter). _q4kRepackSoaKernel, _matvecQ4KSoaKernel, _matvecQ4KN2SoaKernel, _matvecQ4KGemmNSoaKernel, _dequantQ4KF16Kernel, _dequantQ4KF16SoaKernel, - _dequantQ6KF16Kernel, _dequantQ5KF16Kernel, // #162 + _dequantQ6KF16Kernel, _dequantQ6KF16SoaKernel, _dequantQ5KF16Kernel, // #162/#204 _dequantQ40F16Kernel, _headNormPureBatchedKernel, _matvecQ40Dp4aKernel, // #124 _mmqQ40Kernel, _q40RepackSoaKernel, _mmqQ40SoaKernel, _matvecQ40SoaKernel, // #124/#173 _matvecQ40Dp4aSoaKernel, _dequantQ40F16SoaKernel, // #124/#173 @@ -6248,6 +6343,8 @@ private void ForceEagerJit() _addBiasRowsKernel, _mmqQ4kSoaActsN16Kernel, // #201 _mmqQ4kSoaActsN16Bm32Kernel, // #205 + _q6kRepackSoaKernel, // #204 + _mmqQ6kSoaActsN16Kernel, _mmqQ6kSoaActsN16Bm32Kernel, // #204 ]; foreach (nint k in kernels) { @@ -6258,7 +6355,7 @@ private void ForceEagerJit() ReadOnlySpan wsKernelSets = [ _matvecF32WsKernels, _matvecQ4KWsKernels, _matvecQ4KWsSoaKernels, _matvecQ5KWsKernels, _matvecQ6KWsKernels, _matvecQ80WsKernels, _matvecQ80WsSoaKernels, - _matvecQ6KWsSwKernels, + _matvecQ6KWsSwKernels, _matvecQ6KWsSoaKernels, ]; foreach (nint[] set in wsKernelSets) foreach (nint k in set) @@ -6313,8 +6410,10 @@ private void LoadKernelFunctions() _matvecQ4KKernel = GetKernelFunc("llm_matvec_q4k"); _matvecQ4KSoaKernel = GetKernelFunc("llm_matvec_q4k_soa"); _q4kRepackSoaKernel = GetKernelFunc("llm_q4k_repack_soa"); + _q6kRepackSoaKernel = GetKernelFunc("llm_q6k_repack_soa"); // #204 _matvecQ5KKernel = GetKernelFunc("llm_matvec_q5k"); _matvecQ6KKernel = GetKernelFunc("llm_matvec_q6k"); + _matvecQ6KSoaKernel = GetKernelFunc("llm_matvec_q6k_soa"); // #204 _matvecQ40Kernel = GetKernelFunc("llm_matvec_q4_0"); _matvecQ40Dp4aKernel = GetKernelFunc("llm_matvec_q4_0_dp4a"); // #124 _matvecQ80Kernel = GetKernelFunc("llm_matvec_q8_0"); @@ -6324,11 +6423,13 @@ private void LoadKernelFunctions() _matvecQ4KN2SoaKernel = GetKernelFunc("llm_matvec_q4k_n2_soa"); _matvecQ5KN2Kernel = GetKernelFunc("llm_matvec_q5k_n2"); _matvecQ6KN2Kernel = GetKernelFunc("llm_matvec_q6k_n2"); + _matvecQ6KN2SoaKernel = GetKernelFunc("llm_matvec_q6k_n2_soa"); // #204 _matvecF32GemmNKernel = GetKernelFunc("llm_matvec_f32_gemm_n"); _matvecQ4KGemmNKernel = GetKernelFunc("llm_matvec_q4k_gemm_n"); _matvecQ4KGemmNSoaKernel = GetKernelFunc("llm_matvec_q4k_gemm_n_soa"); _matvecQ5KGemmNKernel = GetKernelFunc("llm_matvec_q5k_gemm_n"); _matvecQ6KGemmNKernel = GetKernelFunc("llm_matvec_q6k_gemm_n"); + _matvecQ6KGemmNSoaKernel = GetKernelFunc("llm_matvec_q6k_gemm_n_soa"); // #204 _matvecQ80GemmNKernel = GetKernelFunc("llm_matvec_q8_0_gemm_n"); for (int v = 0; v < CudaWsKernels.Variants.Length; v++) // #194 { @@ -6341,13 +6442,17 @@ private void LoadKernelFunctions() _matvecQ80WsKernels[v] = GetKernelFunc($"llm_matvec_q8_0_ws_n{nt}"); _matvecQ80WsSoaKernels[v] = GetKernelFunc($"llm_matvec_q8_0_ws_soa_n{nt}"); _matvecQ6KWsSwKernels[v] = GetKernelFunc($"llm_matvec_q6k_ws_sw_n{nt}"); // #201 + _matvecQ6KWsSoaKernels[v] = GetKernelFunc($"llm_matvec_q6k_ws_soa_n{nt}"); // #204 } _mmqQ4kSoaActsN16Kernel = GetKernelFunc("llm_mmq_q4k_soa_acts_n16"); // #201 _mmqQ4kSoaActsN16Bm32Kernel = GetKernelFunc("llm_mmq_q4k_soa_acts_n16_bm32"); // #205 + _mmqQ6kSoaActsN16Kernel = GetKernelFunc("llm_mmq_q6k_soa_acts_n16"); // #204 + _mmqQ6kSoaActsN16Bm32Kernel = GetKernelFunc("llm_mmq_q6k_soa_acts_n16_bm32"); // #204 _dequantQ80F16Kernel = GetKernelFunc("llm_dequant_q8_0_to_f16"); _dequantQ4KF16Kernel = GetKernelFunc("llm_dequant_q4k_to_f16"); _dequantQ4KF16SoaKernel = GetKernelFunc("llm_dequant_q4k_to_f16_soa"); _dequantQ6KF16Kernel = GetKernelFunc("llm_dequant_q6k_to_f16"); + _dequantQ6KF16SoaKernel = GetKernelFunc("llm_dequant_q6k_to_f16_soa"); // #204 _dequantQ5KF16Kernel = GetKernelFunc("llm_dequant_q5k_to_f16"); _dequantQ40F16Kernel = GetKernelFunc("llm_dequant_q4_0_to_f16"); // #124 _headNormPureBatchedKernel = GetKernelFunc("llm_head_norm_pure_batched"); // #124 @@ -6835,6 +6940,7 @@ public void Dispose() _devPtrs.Clear(); _soaHandles.Clear(); // #149 _soaQ4kHandles.Clear(); // #156 + _soaQ6kHandles.Clear(); // #204 _soaQ40Handles.Clear(); // #124/#173 _pool.Dispose(); diff --git a/src/SharpInference.Cuda/CudaTextKernels.cs b/src/SharpInference.Cuda/CudaTextKernels.cs index 557e759..1195131 100644 --- a/src/SharpInference.Cuda/CudaTextKernels.cs +++ b/src/SharpInference.Cuda/CudaTextKernels.cs @@ -1631,6 +1631,65 @@ __device__ __forceinline__ void sharpi_q8_append_one( for (int i = 0; i < 128; i++) qDst[i] = src[srcOff + 16 + i]; } +// ── #204 Q6_K → SoA repack for the int8 decode-MMQ tile ───────────────────── +// One-time repack of an interleaved Q6_K weight [rows × nb × 210 B] into a SoA +// layout the llm_mmq_q6k_soa_acts_n16 tile consumes: +// [Q total*256][S total*16][D total*4] (total = rows*nb super-blocks) +// Q stores, for each natural element e of the super-block, the SIGNED int8 weight +// (q6(e) − 32) at byte position e — i.e. the SAME natural-order layout the shared +// Q8_1 activation uses (sub-block sb = e>>5, word (e&31)>>2, byte (e&31)&3), so the +// kernel's a-fragment load mirrors the Q4_K tile's exactly (4 consecutive int8 per +// word). q6 ∈ [0,63] → q6−32 ∈ [−32,31] fits signed int8. S stores the 16 int8 +// scales verbatim (bytes 192..207); D stores {fp16 d (bytes 208..209), 0 pad}. The +// reconstruction mirrors llm_dequant_q6k_to_f16 / DotQ6K exactly, so the int8 weight +// bytes are bit-identical to the matvec's pre-multiply (q − 32). One thread / super-block. +extern ""C"" __global__ void llm_q6k_repack_soa( + const unsigned char* __restrict__ src, // interleaved, 210 B/super-block + unsigned char* __restrict__ dst, // SoA [Q][S][D] + int rows, int cols) +{ + long sb = (long)blockIdx.x * blockDim.x + threadIdx.x; + int nb = cols >> 8; + long total = (long)rows * nb; + if (sb >= total) return; + + long srcOff = sb * 210L; + signed char* qDst = (signed char*)(dst + sb * 256L); + unsigned char* sDst = dst + total * 256L + sb * 16L; + unsigned char* dDst = dst + total * (256L + 16L) + sb * 4L; + + const unsigned char* ql = src + srcOff; // [0:128] + const unsigned char* qh = src + srcOff + 128; // [128:192] + + // 16 int8 scales (192..207) → S region, verbatim. + #pragma unroll + for (int i = 0; i < 16; i++) sDst[i] = src[srcOff + 192 + i]; + // {d (208,209), 0, 0} → D region. + dDst[0] = src[srcOff + 208]; dDst[1] = src[srcOff + 209]; + dDst[2] = 0; dDst[3] = 0; + + // Reconstruct all 256 (q6 − 32) signed int8, element e → byte e (natural order). + // group = e>>5 (0..7), lane = e&31; the ql/qh switch is identical to + // llm_dequant_q6k_to_f16 (same byte/shift per group). + #pragma unroll + for (int e = 0; e < 256; e++) { + int lane = e & 31; + int group = e >> 5; + unsigned int qlb, qhb; int q6; + switch (group) { + case 0: qlb = ql[ 0 + lane]; qhb = qh[ 0 + lane]; q6 = (int)((qlb & 0xFu) | (((qhb >> 0) & 3u) << 4)); break; + case 1: qlb = ql[ 32 + lane]; qhb = qh[ 0 + lane]; q6 = (int)((qlb & 0xFu) | (((qhb >> 2) & 3u) << 4)); break; + case 2: qlb = ql[ 0 + lane]; qhb = qh[ 0 + lane]; q6 = (int)(((qlb >> 4) & 0xFu) | (((qhb >> 4) & 3u) << 4)); break; + case 3: qlb = ql[ 32 + lane]; qhb = qh[ 0 + lane]; q6 = (int)(((qlb >> 4) & 0xFu) | (((qhb >> 6) & 3u) << 4)); break; + case 4: qlb = ql[ 64 + lane]; qhb = qh[32 + lane]; q6 = (int)((qlb & 0xFu) | (((qhb >> 0) & 3u) << 4)); break; + case 5: qlb = ql[ 96 + lane]; qhb = qh[32 + lane]; q6 = (int)((qlb & 0xFu) | (((qhb >> 2) & 3u) << 4)); break; + case 6: qlb = ql[ 64 + lane]; qhb = qh[32 + lane]; q6 = (int)(((qlb >> 4) & 0xFu) | (((qhb >> 4) & 3u) << 4)); break; + default: qlb = ql[ 96 + lane]; qhb = qh[32 + lane]; q6 = (int)(((qlb >> 4) & 0xFu) | (((qhb >> 6) & 3u) << 4)); break; + } + qDst[e] = (signed char)(q6 - 32); + } +} + // ── MatVec Q6_K ──────────────────────────────────────────────────────────── // Q6_K block (210 bytes per 256 elements): // [0:128] ql — lower 4-bit nibbles (two 64-byte halves) @@ -1701,6 +1760,69 @@ __device__ __forceinline__ void sharpi_q8_append_one( if (lane == 0) output[row] = result; } +// ── MatVec Q6_K — SoA (#204) ──────────────────────────────────────────────── +// Bit-identical clone of llm_matvec_q6k over the scale-pre-unpacked SoA layout +// (llm_q6k_repack_soa): [Q total*256 signed-int8 (q6−32)][S total*16 int8 scales] +// [D total*4 {fp16 d, 0, 0}], super-block index g = row*nb + block. Element e of +// the super-block is the signed weight Q[g*256 + e]; its scale is S[g*16 + (e>>4)]; +// d is D[g*4]. The matvec reads the SAME element order (group 0..7, lane) and forms +// the SAME sc·(q6−32)·input products in the SAME reduction order as the AoS kernel, +// so the output is bit-identical to llm_matvec_q6k (no bit-unpacking, no scale gather). +extern ""C"" __global__ void llm_matvec_q6k_soa( + const unsigned char* __restrict__ weights, // SoA [Q][S][D] + const float* __restrict__ input, + float* __restrict__ output, + int rows, int cols) +{ + const int N_ROWS = 8; + const int THREADS_PER_ROW = 32; + unsigned int tid = threadIdx.x; + int row_in_wg = (int)tid / THREADS_PER_ROW; + int lane = (int)tid & (THREADS_PER_ROW - 1); + int row = (int)blockIdx.x * N_ROWS + row_in_wg; + if (row >= rows) return; + + int num_blocks = cols >> 8; + long total_sb = (long)rows * num_blocks; + const signed char* qReg = (const signed char*)weights; // [Q] total*256 B + const signed char* sReg = (const signed char*)weights + total_sb * 256L; // [S] total*16 B + const unsigned char* dReg = (const unsigned char*)weights + total_sb * (256L + 16L); // [D] total*4 B + + long isc = (long)(lane >> 4); // scale half within the 32-element group + + float acc = 0.f; + + for (int block = 0; block < num_blocks; block++) { + long g = (long)row * num_blocks + block; + const signed char* q = qReg + g * 256L; // 256 signed int8 (q6-32) + const signed char* s = sReg + g * 16L; // 16 int8 scales + unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + float d = sharpi_fp16_to_fp32(dbits); + + float sc0 = d * (float)s[ 0 + isc]; + float sc1 = d * (float)s[ 2 + isc]; + float sc2 = d * (float)s[ 4 + isc]; + float sc3 = d * (float)s[ 6 + isc]; + float sc4 = d * (float)s[ 8 + isc]; + float sc5 = d * (float)s[10 + isc]; + float sc6 = d * (float)s[12 + isc]; + float sc7 = d * (float)s[14 + isc]; + + int base_elem = block * 256; + acc += sc0 * (float)q[ lane] * input[base_elem + lane]; + acc += sc1 * (float)q[ 32 + lane] * input[base_elem + 32 + lane]; + acc += sc2 * (float)q[ 64 + lane] * input[base_elem + 64 + lane]; + acc += sc3 * (float)q[ 96 + lane] * input[base_elem + 96 + lane]; + acc += sc4 * (float)q[128 + lane] * input[base_elem + 128 + lane]; + acc += sc5 * (float)q[160 + lane] * input[base_elem + 160 + lane]; + acc += sc6 * (float)q[192 + lane] * input[base_elem + 192 + lane]; + acc += sc7 * (float)q[224 + lane] * input[base_elem + 224 + lane]; + } + + float result = sharpi_warp_reduce_sum(acc); + if (lane == 0) output[row] = result; +} + // ── MatVec Q5_K ──────────────────────────────────────────────────────────── // Q5_K block (176 bytes per 256 elements): // [0:2] fp16 d — super-block scale @@ -4256,6 +4378,43 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( } } +// ── Dequant Q6_K → FP16 SoA (#204) ────────────────────────────────────────── +// Bit-identical clone of llm_dequant_q6k_to_f16 over the SoA layout (see +// llm_matvec_q6k_soa). Thread e owns weight column e of the super-block: scale +// S[g*16 + group*2 + (lane>>4)], weight (q6-32) = Q[g*256 + e], d = D[g*4]. +// val = d·scale·(q6-32), fp16-rounded — same lossy step + same out index as the AoS kernel. +extern ""C"" __global__ void llm_dequant_q6k_to_f16_soa( + const unsigned char* __restrict__ weights, // SoA [Q][S][D] + unsigned short* __restrict__ out, // [rows * cols] fp16 + int rows, int cols) +{ + int row = (int)blockIdx.x; + if (row >= rows) return; + int num_blocks = cols >> 8; // cols / 256 + long out_row = (long)row * (long)cols; + long total_sb = (long)rows * num_blocks; + const signed char* qReg = (const signed char*)weights; + const signed char* sReg = (const signed char*)weights + total_sb * 256L; + const unsigned char* dReg = (const unsigned char*)weights + total_sb * (256L + 16L); + + for (int block = 0; block < num_blocks; block++) { + long g = (long)row * num_blocks + block; + const signed char* q = qReg + g * 256L; + const signed char* s = sReg + g * 16L; + unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + float d = sharpi_fp16_to_fp32(dbits); + + for (int e = (int)threadIdx.x; e < 256; e += (int)blockDim.x) { + int lane = e & 31; + int group = e >> 5; // 0..7 + int isc = lane >> 4; // 0 or 1 + float sc = d * (float)s[group * 2 + isc]; + float val = sc * (float)q[e]; + out[out_row + (long)block * 256 + e] = (unsigned short)sharpi_fp32_to_fp16(val); + } + } +} + // ── Dequant Q5_K → FP16 for cuBLAS prefill GEMM (issue #162) ──────────────── // Same motivation as the Q6_K dequant: Q5_K_M mixes keep q/k/o/gate/up in Q5_K, which // otherwise fell to the per-token GEMM-N matvec. Element decode mirrors llm_matvec_q5k: @@ -4664,6 +4823,89 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( } } +// ── MatVec Q6_K — N=2 SoA (#204) ──────────────────────────────────────────── +// Bit-identical clone of llm_matvec_q6k_n2 over the SoA layout (see +// llm_matvec_q6k_soa). Same per-element dequant + reduction order for both inputs. +extern ""C"" __global__ void llm_matvec_q6k_n2_soa( + const unsigned char* __restrict__ weights, // SoA [Q][S][D] + const float* __restrict__ input_a, + const float* __restrict__ input_b, + float* __restrict__ output_a, + float* __restrict__ output_b, + int rows, int cols) +{ + const int N_ROWS = 8; + const int THREADS_PER_ROW = 32; + unsigned int tid = threadIdx.x; + int row_in_wg = (int)tid / THREADS_PER_ROW; + int lane = (int)tid & (THREADS_PER_ROW - 1); + int row = (int)blockIdx.x * N_ROWS + row_in_wg; + if (row >= rows) return; + + int num_blocks = cols >> 8; + long total_sb = (long)rows * num_blocks; + const signed char* qReg = (const signed char*)weights; + const signed char* sReg = (const signed char*)weights + total_sb * 256L; + const unsigned char* dReg = (const unsigned char*)weights + total_sb * (256L + 16L); + + long isc = (long)(lane >> 4); + + float acc_a = 0.f, acc_b = 0.f; + + for (int block = 0; block < num_blocks; block++) { + long g = (long)row * num_blocks + block; + const signed char* q = qReg + g * 256L; + const signed char* s = sReg + g * 16L; + unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + float d = sharpi_fp16_to_fp32(dbits); + + float sc0 = d * (float)s[ 0 + isc]; + float sc1 = d * (float)s[ 2 + isc]; + float sc2 = d * (float)s[ 4 + isc]; + float sc3 = d * (float)s[ 6 + isc]; + float sc4 = d * (float)s[ 8 + isc]; + float sc5 = d * (float)s[10 + isc]; + float sc6 = d * (float)s[12 + isc]; + float sc7 = d * (float)s[14 + isc]; + + int base_elem = block * 256; + + float q0 = sc0 * (float)q[ lane]; + float q1 = sc1 * (float)q[ 32 + lane]; + float q2 = sc2 * (float)q[ 64 + lane]; + float q3 = sc3 * (float)q[ 96 + lane]; + float q4 = sc4 * (float)q[128 + lane]; + float q5 = sc5 * (float)q[160 + lane]; + float q6 = sc6 * (float)q[192 + lane]; + float q7 = sc7 * (float)q[224 + lane]; + + acc_a += q0 * input_a[base_elem + lane]; + acc_a += q1 * input_a[base_elem + 32 + lane]; + acc_a += q2 * input_a[base_elem + 64 + lane]; + acc_a += q3 * input_a[base_elem + 96 + lane]; + acc_a += q4 * input_a[base_elem + 128 + lane]; + acc_a += q5 * input_a[base_elem + 160 + lane]; + acc_a += q6 * input_a[base_elem + 192 + lane]; + acc_a += q7 * input_a[base_elem + 224 + lane]; + + acc_b += q0 * input_b[base_elem + lane]; + acc_b += q1 * input_b[base_elem + 32 + lane]; + acc_b += q2 * input_b[base_elem + 64 + lane]; + acc_b += q3 * input_b[base_elem + 96 + lane]; + acc_b += q4 * input_b[base_elem + 128 + lane]; + acc_b += q5 * input_b[base_elem + 160 + lane]; + acc_b += q6 * input_b[base_elem + 192 + lane]; + acc_b += q7 * input_b[base_elem + 224 + lane]; + } + + float ra = sharpi_warp_reduce_sum(acc_a); + float rb = sharpi_warp_reduce_sum(acc_b); + if (lane == 0) { + output_a[row] = ra; + output_b[row] = rb; + } +} + // ── MatVec Q5_K — N=2 variant (issue #43) ───────────────────────────────── extern ""C"" __global__ void llm_matvec_q5k_n2( const unsigned int* __restrict__ weights, @@ -5057,6 +5299,67 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( if (lane == 0) output_all[(long)token * (long)rows + row] = result; } +// ── MatVec Q6_K GEMM-N SoA (#204) ─────────────────────────────────────────── +// Bit-identical clone of llm_matvec_q6k_gemm_n over the SoA layout (see +// llm_matvec_q6k_soa). Same per-(row,token) reduction order, so a (rows, n_tok) +// launch is bit-identical to n_tok sequential llm_matvec_q6k_soa calls. +extern ""C"" __global__ void llm_matvec_q6k_gemm_n_soa( + const unsigned char* __restrict__ weights, // SoA [Q][S][D] + const float* __restrict__ input_all, // [n_tok][cols] + float* __restrict__ output_all, // [n_tok][rows] + int rows, int cols, int n_tok) +{ + const int N_ROWS = 8; + const int THREADS_PER_ROW = 32; + unsigned int tid = threadIdx.x; + int row_in_wg = (int)tid / THREADS_PER_ROW; + int lane = (int)tid & (THREADS_PER_ROW - 1); + int row = (int)blockIdx.x * N_ROWS + row_in_wg; + int token = (int)blockIdx.y; + if (row >= rows || token >= n_tok) return; + + const float* input = input_all + (long)token * (long)cols; + int num_blocks = cols >> 8; + long total_sb = (long)rows * num_blocks; + const signed char* qReg = (const signed char*)weights; + const signed char* sReg = (const signed char*)weights + total_sb * 256L; + const unsigned char* dReg = (const unsigned char*)weights + total_sb * (256L + 16L); + + long isc = (long)(lane >> 4); + + float acc = 0.f; + + for (int block = 0; block < num_blocks; block++) { + long g = (long)row * num_blocks + block; + const signed char* q = qReg + g * 256L; + const signed char* s = sReg + g * 16L; + unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + float d = sharpi_fp16_to_fp32(dbits); + + float sc0 = d * (float)s[ 0 + isc]; + float sc1 = d * (float)s[ 2 + isc]; + float sc2 = d * (float)s[ 4 + isc]; + float sc3 = d * (float)s[ 6 + isc]; + float sc4 = d * (float)s[ 8 + isc]; + float sc5 = d * (float)s[10 + isc]; + float sc6 = d * (float)s[12 + isc]; + float sc7 = d * (float)s[14 + isc]; + + int base_elem = block * 256; + acc += sc0 * (float)q[ lane] * input[base_elem + lane]; + acc += sc1 * (float)q[ 32 + lane] * input[base_elem + 32 + lane]; + acc += sc2 * (float)q[ 64 + lane] * input[base_elem + 64 + lane]; + acc += sc3 * (float)q[ 96 + lane] * input[base_elem + 96 + lane]; + acc += sc4 * (float)q[128 + lane] * input[base_elem + 128 + lane]; + acc += sc5 * (float)q[160 + lane] * input[base_elem + 160 + lane]; + acc += sc6 * (float)q[192 + lane] * input[base_elem + 192 + lane]; + acc += sc7 * (float)q[224 + lane] * input[base_elem + 224 + lane]; + } + + float result = sharpi_warp_reduce_sum(acc); + if (lane == 0) output_all[(long)token * (long)rows + row] = result; +} + // ── MatVec Q5_K GEMM-N (issue #119) ──────────────────────────────────────── // Batched (token-dimension) clone of `llm_matvec_q5k`: F32 input, per-element // Q5_K weight decode, 8 rows/block × 32 threads/row warp reduce. Per-lane diff --git a/src/SharpInference.Cuda/CudaWsKernels.cs b/src/SharpInference.Cuda/CudaWsKernels.cs index 9360fef..147f575 100644 --- a/src/SharpInference.Cuda/CudaWsKernels.cs +++ b/src/SharpInference.Cuda/CudaWsKernels.cs @@ -55,24 +55,36 @@ internal static class CudaWsKernels private static string Build() { - var sb = new System.Text.StringBuilder(Template.Length * Variants.Length + DecodeMmqTemplate.Length * 2); + var sb = new System.Text.StringBuilder( + Template.Length * Variants.Length + + (DecodeMmqTemplate.Length + DecodeMmqQ6KTemplate.Length) * 2); foreach (int nt in Variants) sb.Append(Template.Replace("__NT__", nt.ToString(System.Globalization.CultureInfo.InvariantCulture))); // Decode MMQ, two row-tile sizes (#205). BM/16 row-strips: BM=64 → wr = warp & 3, // wc = warp >> 2; BM=32 → wr = warp & 1, wc = warp >> 1. Quant staging is 8 iters // either way (MMQ_BM*32 == 8*threads); act staging is (16*64)/threads = 4 / 8 iters. - sb.Append(EmitDecodeMmq(bm: 64, threads: 256, actJ: 4, wrMask: 3, wcShift: 2, suffix: "")); - sb.Append(EmitDecodeMmq(bm: 32, threads: 128, actJ: 8, wrMask: 1, wcShift: 1, suffix: "_bm32")); + sb.Append(EmitDecodeMmq(DecodeMmqTemplate, bm: 64, threads: 256, actJ: 4, wrMask: 3, wcShift: 2, suffix: "")); + sb.Append(EmitDecodeMmq(DecodeMmqTemplate, bm: 32, threads: 128, actJ: 8, wrMask: 1, wcShift: 1, suffix: "_bm32")); + // #204 Q6_K decode MMQ: same two row-tiles. The Q region is 256 int8 / super-block + // (64 words/row) — quant staging is MMQ_BM*64/threads = 16 / 16 iters. The SoA weight + // is now the ONLY copy of the Q6_K weight (RepackQ6KSoa frees the AoS), so there is no + // AoS-direct decode-MMQ variant — every Q6_K reader is SoA-aware. + sb.Append(EmitDecodeMmq(DecodeMmqQ6KTemplate, bm: 64, threads: 256, actJ: 4, wrMask: 3, wcShift: 2, suffix: "", + wQuantJ: 16)); + sb.Append(EmitDecodeMmq(DecodeMmqQ6KTemplate, bm: 32, threads: 128, actJ: 8, wrMask: 1, wcShift: 1, suffix: "_bm32", + wQuantJ: 16)); return sb.ToString(); } - private static string EmitDecodeMmq(int bm, int threads, int actJ, int wrMask, int wcShift, string suffix) + private static string EmitDecodeMmq(string template, int bm, int threads, int actJ, int wrMask, int wcShift, + string suffix, int wQuantJ = 0) { static string S(int v) => v.ToString(System.Globalization.CultureInfo.InvariantCulture); - return DecodeMmqTemplate + return template .Replace("__BM__", S(bm)) .Replace("__NTHREADS__", S(threads)) .Replace("__ACTJ__", S(actJ)) + .Replace("__WQUANTJ__", S(wQuantJ)) .Replace("__WRMASK__", S(wrMask)) .Replace("__WCSHIFT__", S(wcShift)) .Replace("__SUF__", suffix); @@ -293,6 +305,88 @@ private static string EmitDecodeMmq(int bm, int threads, int actJ, int wrMask, i } } +// Q6_K SoA (#204): bit-identical to llm_matvec_q6k_ws_n* (and to its scale-word +// twin llm_matvec_q6k_ws_sw_n*) over the scale-pre-unpacked SoA layout (see +// llm_matvec_q6k_soa). The 16 int8 scales are already separate, so the scale-word +// trick is moot — read S[g*16 + 2k + (lane>>4)] directly; the (q6−32) int8 weight is +// Q[g*256 + e], d is D[g*4]. Same 8-element-group reduction order + same token loop, +// so the output is bit-identical to both AoS Q6_K WS kernels. +extern ""C"" __global__ void llm_matvec_q6k_ws_soa_n__NT__( + const unsigned char* __restrict__ weights, // SoA [Q][S][D] + const float* __restrict__ input_all, // [n_tok][cols] + float* __restrict__ output_all, // [n_tok][rows] + int rows, int cols, int n_tok) +{ + const int NT = __NT__; + const int N_ROWS = 8; + const int THREADS_PER_ROW = 32; + unsigned int tid = threadIdx.x; + int row_in_wg = (int)tid / THREADS_PER_ROW; + int lane = (int)tid & (THREADS_PER_ROW - 1); + int row = (int)blockIdx.x * N_ROWS + row_in_wg; + if (row >= rows) return; + + int num_blocks = cols >> 8; + long total_sb = (long)rows * num_blocks; + const signed char* qReg = (const signed char*)weights; + const signed char* sReg = (const signed char*)weights + total_sb * 256L; + const unsigned char* dReg = (const unsigned char*)weights + total_sb * (256L + 16L); + + long isc = (long)(lane >> 4); + + float acc[NT]; + #pragma unroll + for (int t = 0; t < NT; t++) acc[t] = 0.f; + + for (int block = 0; block < num_blocks; block++) { + long g = (long)row * num_blocks + block; + const signed char* q = qReg + g * 256L; + const signed char* s = sReg + g * 16L; + unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + float d = sharpi_fp16_to_fp32(dbits); + + float sc0 = d * (float)s[ 0 + isc]; + float sc1 = d * (float)s[ 2 + isc]; + float sc2 = d * (float)s[ 4 + isc]; + float sc3 = d * (float)s[ 6 + isc]; + float sc4 = d * (float)s[ 8 + isc]; + float sc5 = d * (float)s[10 + isc]; + float sc6 = d * (float)s[12 + isc]; + float sc7 = d * (float)s[14 + isc]; + + float w0 = sc0 * (float)q[ lane]; + float w1 = sc1 * (float)q[ 32 + lane]; + float w2 = sc2 * (float)q[ 64 + lane]; + float w3 = sc3 * (float)q[ 96 + lane]; + float w4 = sc4 * (float)q[128 + lane]; + float w5 = sc5 * (float)q[160 + lane]; + float w6 = sc6 * (float)q[192 + lane]; + float w7 = sc7 * (float)q[224 + lane]; + + int base_elem = block * 256; + #pragma unroll + for (int t = 0; t < NT; t++) + if (t < n_tok) { + const float* input = input_all + (long)t * (long)cols; + acc[t] += w0 * input[base_elem + lane]; + acc[t] += w1 * input[base_elem + 32 + lane]; + acc[t] += w2 * input[base_elem + 64 + lane]; + acc[t] += w3 * input[base_elem + 96 + lane]; + acc[t] += w4 * input[base_elem + 128 + lane]; + acc[t] += w5 * input[base_elem + 160 + lane]; + acc[t] += w6 * input[base_elem + 192 + lane]; + acc[t] += w7 * input[base_elem + 224 + lane]; + } + } + + #pragma unroll + for (int t = 0; t < NT; t++) + if (t < n_tok) { + float result = sharpi_warp_reduce_sum(acc[t]); + if (lane == 0) output_all[(long)t * (long)rows + row] = result; + } +} + // Q5_K: hoists the token-invariant (d1*q - dm1) per element pair — the same // fma llm_matvec_q5k_gemm_n contracts to before the activation multiply. extern ""C"" __global__ void llm_matvec_q5k_ws_n__NT__( @@ -916,5 +1010,175 @@ private static string EmitDecodeMmq(int bm, int threads, int actJ, int wrMask, i } #undef MMQ_BM #undef MMQ_BN +"; + + /// + /// #204 decode MMQ for Q6_K weights: the int8 tensor-core analogue of + /// for the Q6_K trunk shapes (Qwen3-8B Q4_K_M keeps + /// half of ffn_down + attn_v + lm-head in Q6_K), which otherwise fall to the bit-exact + /// llm_matvec_q6k_ws_sw matvec (~2.7-3× the weight-streaming floor at N=8). + /// + /// The SoA weight (llm_q6k_repack_soa) stores the signed int8 + /// (q6 − 32) for natural element e at byte e — the SAME natural order the + /// shared Q8_1 activation uses — so the a-fragment load and the b-fragment (activation) + /// load mirror the Q4_K tile exactly. Q6_K is SYMMETRIC: no min term, the −32 offset + /// is baked into the int8 weight. The only structural change vs the Q4_K tile is the + /// per-16-element scale boundary inside each 32-element sub-block: registers a0/a1 of + /// the m16n8k32 fragment hold k∈[0,16), a2/a3 hold k∈[16,32) — a clean split at that + /// boundary — so each sub-block runs TWO masked mmas (mma1 = {a0,a1,0,0}, mma2 = + /// {0,0,a2,a3}) and accumulates c1·dwG0·dC + c2·dwG1·dC with dwG0 = d·scales[2·sb], + /// dwG1 = d·scales[2·sb+1]. Argmax-stable (both operands int8-quantized), not bit-exact — + /// same contract as the Q4_K decode/prefill MMQ. + /// + private const string DecodeMmqQ6KTemplate = @" +// ── #204 decode MMQ: Q6_K SoA weights × SoA Q8_1 activations, BN=16 ────────── +// grid = (ceil(rows/MMQ_BM), ceil(n_tok/16)), block = __NTHREADS__. K-step is one +// 256-element super-block: stage the raw tile (MMQ_BM×256 int8 quant words + 16-B +// scales + {d} + 16×64 act words + 16×8 act {d,s}) with coalesced copies, then 8 +// sub-blocks × TWO m16n8k32 s8 mma per warp between one barrier pair. The two +// half-fragment mmas split the per-16-element Q6_K scale boundary; the −32 offset is +// in the int8 weight so there is no min term. Argmax-stable (both operands int8). +#define MMQ_BM __BM__ +#define MMQ_BN 16 +extern ""C"" __global__ void __launch_bounds__(__NTHREADS__) llm_mmq_q6k_soa_acts_n16__SUF__( + const unsigned int* __restrict__ weights, // SoA [Q total*256][S total*16][D total*4] + const unsigned int* __restrict__ y_qs, // SoA Q8_1 quants [n_tok × sub_total × 32 B] + const unsigned int* __restrict__ y_ds, // SoA Q8_1 scales [n_tok × sub_total] {d,s} + float* __restrict__ output, // [n_tok × rows] fp32 + int rows, int cols, int n_tok) +{ + // 64 quant words/row (256 int8). +1 word stride avoids the 8-way bank conflict an + // unpadded 64-word stride would land on the 8 grp-lanes of a fragment read. + __shared__ unsigned int sWq[MMQ_BM * 65]; // raw int8 quant words [row][64+1] + __shared__ unsigned int sWs[MMQ_BM * 4]; // 16 int8 scales [row][4] + __shared__ unsigned int sWd[MMQ_BM]; // {d, 0} [row] + __shared__ unsigned int sYq[MMQ_BN * 65]; // raw act words [tok][64+1] + __shared__ unsigned int sYd[MMQ_BN * 8]; // raw act {d,s} [tok][8] + + int row_block = (int)blockIdx.x * MMQ_BM; + int tok_block = (int)blockIdx.y * MMQ_BN; + int nb_super = cols >> 8; + int sub_total = cols >> 5; + long total_sb = (long)rows * nb_super; + + const unsigned int* qReg = weights; // [Q] total*256 B + const unsigned char* sReg = (const unsigned char*)weights + total_sb * 256L; // [S] total*16 B + const unsigned int* dReg = (const unsigned int*)(sReg + total_sb * 16L); // [D] total*4 B + + int tid = (int)threadIdx.x; + int warp = tid >> 5; + int lane = tid & 31; + int grp = lane >> 2; + int tig = lane & 3; + int wr = warp & __WRMASK__; + int wc = warp >> __WCSHIFT__; + int mrow0 = wr * 16; + int ncol0 = wc * 8; + + float acc0 = 0.f, acc1 = 0.f, acc2 = 0.f, acc3 = 0.f; + + int rowA = row_block + mrow0 + grp; + int rowB = rowA + 8; + int tokC0 = tok_block + ncol0 + tig * 2; + int tokC1 = tokC0 + 1; + + for (int ksb = 0; ksb < nb_super; ksb++) { + // Stage the raw super-block tile, fully coalesced. + #pragma unroll + for (int j = 0; j < __WQUANTJ__; j++) { // MMQ_BM rows × 64 quant words + int k = tid + j * __NTHREADS__; + int r = row_block + (k >> 6); + sWq[(k >> 6) * 65 + (k & 63)] = + (r < rows) ? qReg[((long)r * nb_super + ksb) * 64L + (k & 63)] : 0u; + } + { // MMQ_BM rows × 4 scale words (16 int8) + int r = row_block + (tid >> 2); + sWs[tid] = (r < rows) + ? reinterpret_cast(sReg)[((long)r * nb_super + ksb) * 4L + (tid & 3)] + : 0u; + } + if (tid < MMQ_BM) { // MMQ_BM rows × {d,0} + int r = row_block + tid; + sWd[tid] = (r < rows) ? dReg[(long)r * nb_super + ksb] : 0u; + } + #pragma unroll + for (int j = 0; j < __ACTJ__; j++) { // 16 tokens × 64 act words (256 B) + int k = tid + j * __NTHREADS__; + int t = tok_block + (k >> 6); + sYq[(k >> 6) * 65 + (k & 63)] = + (t < n_tok) ? y_qs[((long)t * sub_total + ksb * 8) * 8L + (k & 63)] : 0u; + } + if (tid < MMQ_BN * 8) { // 16 tokens × 8 act {d,s} + int t = tok_block + (tid >> 3); + sYd[tid] = (t < n_tok) ? y_ds[(long)t * sub_total + ksb * 8 + (tid & 7)] : 0u; + } + __syncthreads(); + + const signed char* sWsB = reinterpret_cast(sWs); + #pragma unroll + for (int sb = 0; sb < 8; sb++) { + // a-fragment: 4 consecutive int8 per word, natural-order (mirrors the Q4_K + // tile's word load but reads int8 directly). a0/a1 = k∈[0,16) (words tig of + // the lower 16 lanes), a2/a3 = k∈[16,32) (words tig+4). The 65-word row stride + // is uint; sub-block sb occupies int8 bytes [sb*32 .. sb*32+31] = uint words + // [sb*8 .. sb*8+7] within the row. + int wbaseA = (mrow0 + grp) * 65 + sb * 8; + int wbaseB = (mrow0 + grp + 8) * 65 + sb * 8; + int a0 = (int)sWq[wbaseA + tig]; + int a1 = (int)sWq[wbaseB + tig]; + int a2 = (int)sWq[wbaseA + tig + 4]; + int a3 = (int)sWq[wbaseB + tig + 4]; + + // Q6_K symmetric scales: sub-block sb's first 16 elems use scales[2*sb], + // second 16 use scales[2*sb+1]. d in low 16 bits of sWd. + float dA = sharpi_fp16_to_fp32(sWd[mrow0 + grp] & 0xffffu); + float dB = sharpi_fp16_to_fp32(sWd[mrow0 + grp + 8] & 0xffffu); + float dwA0 = dA * (float)sWsB[(mrow0 + grp) * 16 + 2 * sb]; + float dwA1 = dA * (float)sWsB[(mrow0 + grp) * 16 + 2 * sb + 1]; + float dwB0 = dB * (float)sWsB[(mrow0 + grp + 8) * 16 + 2 * sb]; + float dwB1 = dB * (float)sWsB[(mrow0 + grp + 8) * 16 + 2 * sb + 1]; + + int b0 = (int)sYq[(ncol0 + grp) * 65 + sb * 8 + tig]; + int b1 = (int)sYq[(ncol0 + grp) * 65 + sb * 8 + tig + 4]; + unsigned int dy0 = sYd[(ncol0 + tig * 2) * 8 + sb]; + unsigned int dy1 = sYd[(ncol0 + tig * 2 + 1) * 8 + sb]; + float dC0 = sharpi_fp16_to_fp32(dy0 & 0xffffu); + float dC1 = sharpi_fp16_to_fp32(dy1 & 0xffffu); + + // mma1 = {a0,a1,0,0} × {b0,b1}: Σ_{k=0..15} (q6−32)·q8 (scale group 2·sb). + int zero = 0; + int p0 = 0, p1 = 0, p2 = 0, p3 = 0; + asm( + ""mma.sync.aligned.m16n8k32.row.col.s32.s8.s8.s32 "" + ""{%0,%1,%2,%3}, {%4,%5,%6,%7}, {%8,%9}, {%0,%1,%2,%3};"" + : ""+r""(p0), ""+r""(p1), ""+r""(p2), ""+r""(p3) + : ""r""(a0), ""r""(a1), ""r""(zero), ""r""(zero), ""r""(b0), ""r""(b1)); + // mma2 = {0,0,a2,a3} × {b0,b1}: Σ_{k=16..31} (q6−32)·q8 (scale group 2·sb+1). + int q0 = 0, q1 = 0, q2 = 0, q3 = 0; + asm( + ""mma.sync.aligned.m16n8k32.row.col.s32.s8.s8.s32 "" + ""{%0,%1,%2,%3}, {%4,%5,%6,%7}, {%8,%9}, {%0,%1,%2,%3};"" + : ""+r""(q0), ""+r""(q1), ""+r""(q2), ""+r""(q3) + : ""r""(zero), ""r""(zero), ""r""(a2), ""r""(a3), ""r""(b0), ""r""(b1)); + + acc0 += ((float)p0 * dwA0 + (float)q0 * dwA1) * dC0; + acc1 += ((float)p1 * dwA0 + (float)q1 * dwA1) * dC1; + acc2 += ((float)p2 * dwB0 + (float)q2 * dwB1) * dC0; + acc3 += ((float)p3 * dwB0 + (float)q3 * dwB1) * dC1; + } + __syncthreads(); + } + + if (rowA < rows) { + if (tokC0 < n_tok) output[(long)tokC0 * rows + rowA] = acc0; + if (tokC1 < n_tok) output[(long)tokC1 * rows + rowA] = acc1; + } + if (rowB < rows) { + if (tokC0 < n_tok) output[(long)tokC0 * rows + rowB] = acc2; + if (tokC1 < n_tok) output[(long)tokC1 * rows + rowB] = acc3; + } +} +#undef MMQ_BM +#undef MMQ_BN "; } diff --git a/src/SharpInference.Engine/CudaForwardPass.cs b/src/SharpInference.Engine/CudaForwardPass.cs index aaeb37b..bdfc206 100644 --- a/src/SharpInference.Engine/CudaForwardPass.cs +++ b/src/SharpInference.Engine/CudaForwardPass.cs @@ -391,6 +391,7 @@ public static DType ResolveConfiguredKvDType() => private bool _forceFlashTc1; // #147 A/B: pin the single-warp TC kernel private readonly bool _mmqSoa; // #149: repack 2-D Q8_0 weights to SoA at upload private readonly bool _q4kSoa; // #156: repack 2-D Q4_K weights to scale-unpacked SoA + private readonly bool _q6kSoa; // #204: repack 2-D Q6_K weights to SoA (frees AoS) private int _bpCapacity; // current N the scratch is sized for (0 = none) private Tensor? _bpHidden, _bpResidual, _bpNorm; // [N × embDim] private Tensor? _bpQ, _bpAttnOut; // [N × numHeads*maxHeadDim] @@ -611,6 +612,15 @@ public CudaForwardPass(GgufModel model, CudaBackend gpu, ModelHyperparams hp, // reverts). Dense-only: the MoE Q4_K readers are not SoA-converted, so the // repack is gated on !_isMoE at upload. _q4kSoa = Environment.GetEnvironmentVariable("SHARPI_Q4K_SOA") != "0"; + // Issue #204: repack ALL 2-D Q6_K trunk weights (ffn_down-half / attn_v / lm-head / + // output.weight in Qwen3-8B Q4_K_M) into the SoA [(q6−32) int8][scales][d] layout AND + // free the interleaved AoS copy — the SoA buffer is the only copy (net ~+0.4 GB over AoS, + // no VRAM-doubling companion). Every Q6_K reader is SoA-aware; the big shapes (rows ≥ + // 2048) route to the int8 decode-MMQ tile, attn_v (rows=1024) to the SoA WS matvec. + // Default on; SHARPI_Q6K_SOA=0 keeps the AoS weight everywhere and disables the decode + // MMQ. Dense-only: the MoE Q6_K readers are not SoA-converted, so the repack is gated on + // !_isMoE at upload. + _q6kSoa = Environment.GetEnvironmentVariable("SHARPI_Q6K_SOA") != "0"; _tqEnabled = enableTurboQuant; _tqBits = enableTurboQuant ? tqBits : 0; @@ -4189,6 +4199,24 @@ private Tensor UploadWeight(string name) int rows = (int)info.Dimensions[1]; result = _gpu.RepackQ4KSoa(result, rows, cols); } + // #204: repack ALL 2-D Q6_K trunk weights (Qwen3-8B Q4_K_M keeps half of ffn_down + + // attn_v + lm-head/output.weight in Q6_K) into the scale-pre-unpacked SoA layout and + // FREE the interleaved AoS copy — the SoA buffer becomes the only copy (net ~+0.4 GB + // over AoS, no VRAM-doubling companion). Every Q6_K reader is SoA-aware, so this fits + // the 12 GB card at N=8 (no spill) while the decode-MMQ holds its win. The token + // embedding table is uploaded separately (UploadRaw, read by llm_embed_lookup_q6k — + // not a matmul) so it never reaches here, but guard the name defensively. cols % 256 + // is required (every Q6_K hidden dim satisfies it); the rows ≥ 2048 floor is dropped + // for the REPACK (attn_v rows=1024 is repacked too, just routes to the SoA WS reader — + // the decode-MMQ rows ≥ 2048 eligibility lives in MatMulBatchedDecodeMmq). Dense-only. + else if (_q6kSoa && !_isMoE && info.DType == DType.Q6_K && info.NDimensions == 2 + && !name.Contains("token_embd", StringComparison.Ordinal)) + { + int cols = (int)info.Dimensions[0]; + int rows = (int)info.Dimensions[1]; + if ((cols & 0xff) == 0) + result = _gpu.RepackQ6KSoa(result, rows, cols); + } // #124/#173: same funnelshift-killing SoA repack for 2-D Q4_0 trunk weights // (Gemma 4 12B QAT). cols % 32 required — every Q4_0 hidden dim satisfies it. // Gated on the same SHARPI_MMQ_SOA flag as Q8_0 (#149); dense-only. diff --git a/tests/SharpInference.Tests.ForwardPass/CudaDecodeMmqTests.cs b/tests/SharpInference.Tests.ForwardPass/CudaDecodeMmqTests.cs index 8f6d439..5f2683e 100644 --- a/tests/SharpInference.Tests.ForwardPass/CudaDecodeMmqTests.cs +++ b/tests/SharpInference.Tests.ForwardPass/CudaDecodeMmqTests.cs @@ -120,6 +120,97 @@ public void DecodeMmq_Q4K_Soa_TracksCpuReference() } } + // ── #204 Q6_K decode MMQ ──────────────────────────────────────────────── + private static byte[] BuildQ6KMatrix(int rows, int cols, Random rng) + { + int blocksPerRow = cols / 256; + int bytesPerRow = blocksPerRow * 210; + var bytes = new byte[rows * bytesPerRow]; + for (int r = 0; r < rows; r++) + for (int b = 0; b < blocksPerRow; b++) + { + int off = r * bytesPerRow + b * 210; + // ql[0:128], qh[128:192], scales[192:208] (int8), d (208:210 fp16). + for (int i = 0; i < 192; i++) + bytes[off + i] = (byte)rng.Next(256); + for (int i = 0; i < 16; i++) + bytes[off + 192 + i] = (byte)(sbyte)(rng.Next(127) - 63); // signed int8 scale + float d = (float)(rng.NextDouble() * 0.04 + 0.005); + ushort dHalf = HalfToUshort((Half)d); + bytes[off + 208] = (byte)(dHalf & 0xFF); + bytes[off + 209] = (byte)(dHalf >> 8); + } + return bytes; + } + + private static void RunCaseQ6K(CudaBackend gpu, Tensor gpuW, byte[] weightBytes, + int rows, int cols, int nTok, Random rng) + { + var acts = new float[nTok * cols]; + for (int i = 0; i < acts.Length; i++) + acts[i] = (float)(rng.NextDouble() * 2 - 1); + + int bytesPerRow = (cols / 256) * 210; + var cpuOut = new float[nTok * rows]; + fixed (byte* wPtr = weightBytes) + fixed (float* aPtr = acts) + { + for (int t = 0; t < nTok; t++) + for (int r = 0; r < rows; r++) + cpuOut[t * rows + r] = SimdKernels.DotQ6K(wPtr + r * bytesPerRow, aPtr + t * cols, cols); + } + + var gpuX = gpu.Upload(acts, TensorShape.D1(acts.Length)); + var gpuY = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + + gpu.MatMulBatchedDecodeMmq(gpuY, gpuW, gpuX, nTok, DType.Q6_K); + gpu.Synchronize(); + + var gpuOut = new float[nTok * rows]; + gpu.Download(gpuY, gpuOut); + gpu.Free(gpuX); + gpu.Free(gpuY); + + double sumSq = 0; + for (int i = 0; i < cpuOut.Length; i++) sumSq += (double)cpuOut[i] * cpuOut[i]; + float refRms = (float)Math.Sqrt(sumSq / cpuOut.Length); + + int mismatches = 0; + float maxAbs = 0; + for (int i = 0; i < cpuOut.Length; i++) + { + float diff = MathF.Abs(gpuOut[i] - cpuOut[i]); + maxAbs = MathF.Max(maxAbs, diff); + if (diff > 0.03f * refRms) mismatches++; + } + Assert.True(mismatches <= cpuOut.Length / 100 + 1, + $"Q6_K decode MMQ rows={rows} cols={cols} nTok={nTok} drifted from fp32 reference: " + + $"{mismatches}/{cpuOut.Length} beyond 3% of RMS ({refRms:E3}), maxAbs={maxAbs:E3}."); + } + + [Fact] + public void DecodeMmq_Q6K_Soa_TracksCpuReference() + { + using var gpu = TryCreate(); + if (gpu is null) return; + + // rows ≥ 2048 (eligible). Same #205 dual-tile coverage as the Q4_K case: 2048 (BM=32, + // no tail), 2096 (BM=32 with a non-multiple-of-32 row tail), 8192 (BM=64). Batch sizes + // cover partial token tiles (2, 5), the full tile (16), and the grid.y = 2 round-up (17). + foreach ((int rows, int cols) in new[] { (2048, 512), (2096, 256), (8192, 512) }) + { + var rng = new Random(20260616 + rows * 13 + cols * 5); + byte[] weightBytes = BuildQ6KMatrix(rows, cols, rng); + var gpuWAos = gpu.UploadRaw(weightBytes, TensorShape.D1(weightBytes.Length), DType.Q6_K); + // #204: RepackQ6KSoa now FREES the AoS weight and returns the SoA buffer (the only + // copy); MatMulBatchedDecodeMmq reads it directly via llm_mmq_q6k_soa_acts_n16. + var gpuW = gpu.RepackQ6KSoa(gpuWAos, rows, cols); + foreach (int nTok in new[] { 2, 5, 8, 16, 17 }) + RunCaseQ6K(gpu, gpuW, weightBytes, rows, cols, nTok, rng); + gpu.Free(gpuW); + } + } + /// Ineligible shapes (rows < 2048 here) must take the weight-stationary /// fallback — same fp32-tracking contract holds trivially (the WS path is bit-exact /// to the per-token matvec). diff --git a/tests/SharpInference.Tests.ForwardPass/CudaGemmQ6KTests.cs b/tests/SharpInference.Tests.ForwardPass/CudaGemmQ6KTests.cs index ec1cc6e..d70a689 100644 --- a/tests/SharpInference.Tests.ForwardPass/CudaGemmQ6KTests.cs +++ b/tests/SharpInference.Tests.ForwardPass/CudaGemmQ6KTests.cs @@ -115,4 +115,117 @@ public void MatMulBatchedGemm_Q6K_TracksCpuReference() $"Q6_K GEMM drifted from fp32 reference: {mismatches}/{cpuOut.Length} beyond 4% of RMS ({refRms:E3}), maxAbs={maxAbs:E3}."); } } + + /// #204: every Q6_K reader ported to the SoA layout + /// () must be BIT-IDENTICAL to its AoS counterpart — + /// same value, same reduction order. This test runs each reader (single-token + /// MatMul, two-input MatMulN2, batched GEMM-N MatMulBatched, and the + /// weight-stationary MatMulBatchedWeightStationary) against the SAME weight bytes, + /// once via the AoS-uploaded weight and once via the SoA-repacked weight, and asserts the + /// two GPU outputs are byte-for-byte equal. (The decode-MMQ tile is only argmax-stable, so + /// it is covered separately by against the fp32 reference.) + /// + [Fact] + public void Q6KSoaReaders_AreBitIdenticalToAos() + { + using var gpu = TryCreate(); + if (gpu is null) return; + + // Shapes: attn_v-like (1024 rows → WS/decode-ineligible), ffn_down-half-like (2048), + // a partial-tile (300), and a wider multi-superblock (512 rows × 1024 cols). + foreach ((int rows, int cols) in new[] { (1024, 512), (2048, 256), (300, 512), (512, 1024) }) + { + var rng = new Random(20260616 + rows * 17 + cols * 3); + byte[] weightBytes = BuildQ6KMatrix(rows, cols, rng); + + // Two GPU copies of the SAME bytes: AoS (interleaved) and SoA (repack frees its src). + var aosW = gpu.UploadRaw(weightBytes, TensorShape.D1(weightBytes.Length), DType.Q6_K); + var soaSrc = gpu.UploadRaw(weightBytes, TensorShape.D1(weightBytes.Length), DType.Q6_K); + var soaW = gpu.RepackQ6KSoa(soaSrc, rows, cols); // frees soaSrc, returns the SoA buffer + + // ── Single-token MatMul ── + { + var x = new float[cols]; + for (int i = 0; i < cols; i++) x[i] = (float)(rng.NextDouble() * 2 - 1); + var gx = gpu.Upload(x, TensorShape.D1(cols)); + var aosY = gpu.Allocate(TensorShape.D1(rows)); + var soaY = gpu.Allocate(TensorShape.D1(rows)); + gpu.MatMul(aosY, aosW, gx, DType.Q6_K); + gpu.MatMul(soaY, soaW, gx, DType.Q6_K); + gpu.Synchronize(); + AssertBitIdentical(gpu, aosY, soaY, rows, $"MatMul rows={rows} cols={cols}"); + gpu.Free(gx); gpu.Free(aosY); gpu.Free(soaY); + } + + // ── Two-input MatMulN2 ── + { + var xa = new float[cols]; + var xb = new float[cols]; + for (int i = 0; i < cols; i++) { xa[i] = (float)(rng.NextDouble() * 2 - 1); xb[i] = (float)(rng.NextDouble() * 2 - 1); } + var gxa = gpu.Upload(xa, TensorShape.D1(cols)); + var gxb = gpu.Upload(xb, TensorShape.D1(cols)); + var aosYa = gpu.Allocate(TensorShape.D1(rows)); + var aosYb = gpu.Allocate(TensorShape.D1(rows)); + var soaYa = gpu.Allocate(TensorShape.D1(rows)); + var soaYb = gpu.Allocate(TensorShape.D1(rows)); + gpu.MatMulN2(aosYa, aosYb, aosW, gxa, gxb, DType.Q6_K); + gpu.MatMulN2(soaYa, soaYb, soaW, gxa, gxb, DType.Q6_K); + gpu.Synchronize(); + AssertBitIdentical(gpu, aosYa, soaYa, rows, $"MatMulN2.A rows={rows} cols={cols}"); + AssertBitIdentical(gpu, aosYb, soaYb, rows, $"MatMulN2.B rows={rows} cols={cols}"); + gpu.Free(gxa); gpu.Free(gxb); + gpu.Free(aosYa); gpu.Free(aosYb); gpu.Free(soaYa); gpu.Free(soaYb); + } + + // ── Batched GEMM-N (MatMulBatched) and weight-stationary, across capacities ── + foreach (int nTok in new[] { 1, 2, 4, 8, 16 }) + { + var acts = new float[(long)nTok * cols]; + for (int i = 0; i < acts.Length; i++) acts[i] = (float)(rng.NextDouble() * 2 - 1); + var gx = gpu.Upload(acts, TensorShape.D1(acts.Length)); + + var aosY = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + var soaY = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + gpu.MatMulBatched(aosY, aosW, gx, nTok, DType.Q6_K); + gpu.MatMulBatched(soaY, soaW, gx, nTok, DType.Q6_K); + gpu.Synchronize(); + AssertBitIdentical(gpu, aosY, soaY, nTok * rows, $"MatMulBatched n={nTok} rows={rows} cols={cols}"); + + // WS path (nTok ≥ 2 uses the WS kernels; nTok == 1 delegates to GEMM-N). + var aosWs = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + var soaWs = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + gpu.MatMulBatchedWeightStationary(aosWs, aosW, gx, nTok, DType.Q6_K); + gpu.MatMulBatchedWeightStationary(soaWs, soaW, gx, nTok, DType.Q6_K); + gpu.Synchronize(); + AssertBitIdentical(gpu, aosWs, soaWs, nTok * rows, $"MatMulBatchedWS n={nTok} rows={rows} cols={cols}"); + + // Prefill dequant→fp16→cuBLAS GEMM path (llm_dequant_q6k_to_f16{,_soa}): the SoA + // dequant emits the SAME fp16 weight bytes as the AoS dequant, so the GEMM output + // is bit-identical too (the hot prefill path — covers the SoA dequant kernel). + var aosG = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + var soaG = gpu.Allocate(TensorShape.D1((long)nTok * rows)); + gpu.MatMulBatchedGemm(aosG, aosW, gx, nTok, DType.Q6_K); + gpu.MatMulBatchedGemm(soaG, soaW, gx, nTok, DType.Q6_K); + gpu.Synchronize(); + AssertBitIdentical(gpu, aosG, soaG, nTok * rows, $"MatMulBatchedGemm n={nTok} rows={rows} cols={cols}"); + + gpu.Free(gx); gpu.Free(aosY); gpu.Free(soaY); gpu.Free(aosWs); gpu.Free(soaWs); + gpu.Free(aosG); gpu.Free(soaG); + } + + gpu.Free(aosW); + gpu.Free(soaW); + } + } + + private static void AssertBitIdentical(CudaBackend gpu, Tensor a, Tensor b, int n, string ctx) + { + var ha = new float[n]; + var hb = new float[n]; + gpu.Download(a, ha); + gpu.Download(b, hb); + for (int i = 0; i < n; i++) + Assert.True(BitConverter.SingleToInt32Bits(ha[i]) == BitConverter.SingleToInt32Bits(hb[i]), + $"{ctx}: SoA Q6_K reader not bit-identical to AoS at [{i}]: AoS={ha[i]:R} SoA={hb[i]:R}."); + } } From 21fa24d925e23af6743e15b4f4f1e5ceaf8bfa74 Mon Sep 17 00:00:00 2001 From: Pekka Heikura Date: Tue, 16 Jun 2026 14:37:32 +0300 Subject: [PATCH 2/2] review: single 16-bit d load in Q6_K SoA readers (gemini) dReg is 16-B aligned and g*4 is 4-aligned, so the per-super-block d (fp16) load is a single aligned u16 read instead of two byte loads + shift/or. Bit-identical (same fp16 value); CudaGemmQ6KTests bit-identity + DecodeMmq_Q6K + 48-token SpecDecode greedy parity still green. Co-Authored-By: Claude Opus 4.8 --- src/SharpInference.Cuda/CudaTextKernels.cs | 8 ++++---- src/SharpInference.Cuda/CudaWsKernels.cs | 2 +- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/SharpInference.Cuda/CudaTextKernels.cs b/src/SharpInference.Cuda/CudaTextKernels.cs index 1195131..fe89c99 100644 --- a/src/SharpInference.Cuda/CudaTextKernels.cs +++ b/src/SharpInference.Cuda/CudaTextKernels.cs @@ -1796,7 +1796,7 @@ __device__ __forceinline__ void sharpi_q8_append_one( long g = (long)row * num_blocks + block; const signed char* q = qReg + g * 256L; // 256 signed int8 (q6-32) const signed char* s = sReg + g * 16L; // 16 int8 scales - unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + unsigned int dbits = (unsigned int)(*(const unsigned short*)(dReg + g * 4L)); // #204 review: dReg 16-B aligned, g*4 aligned → single 16-bit d load float d = sharpi_fp16_to_fp32(dbits); float sc0 = d * (float)s[ 0 + isc]; @@ -4401,7 +4401,7 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( long g = (long)row * num_blocks + block; const signed char* q = qReg + g * 256L; const signed char* s = sReg + g * 16L; - unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + unsigned int dbits = (unsigned int)(*(const unsigned short*)(dReg + g * 4L)); // #204 review: dReg 16-B aligned, g*4 aligned → single 16-bit d load float d = sharpi_fp16_to_fp32(dbits); for (int e = (int)threadIdx.x; e < 256; e += (int)blockDim.x) { @@ -4856,7 +4856,7 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( long g = (long)row * num_blocks + block; const signed char* q = qReg + g * 256L; const signed char* s = sReg + g * 16L; - unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + unsigned int dbits = (unsigned int)(*(const unsigned short*)(dReg + g * 4L)); // #204 review: dReg 16-B aligned, g*4 aligned → single 16-bit d load float d = sharpi_fp16_to_fp32(dbits); float sc0 = d * (float)s[ 0 + isc]; @@ -5333,7 +5333,7 @@ __device__ __forceinline__ void sharpi_q4k_scale_min( long g = (long)row * num_blocks + block; const signed char* q = qReg + g * 256L; const signed char* s = sReg + g * 16L; - unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + unsigned int dbits = (unsigned int)(*(const unsigned short*)(dReg + g * 4L)); // #204 review: dReg 16-B aligned, g*4 aligned → single 16-bit d load float d = sharpi_fp16_to_fp32(dbits); float sc0 = d * (float)s[ 0 + isc]; diff --git a/src/SharpInference.Cuda/CudaWsKernels.cs b/src/SharpInference.Cuda/CudaWsKernels.cs index 147f575..3a9dd00 100644 --- a/src/SharpInference.Cuda/CudaWsKernels.cs +++ b/src/SharpInference.Cuda/CudaWsKernels.cs @@ -342,7 +342,7 @@ private static string EmitDecodeMmq(string template, int bm, int threads, int ac long g = (long)row * num_blocks + block; const signed char* q = qReg + g * 256L; const signed char* s = sReg + g * 16L; - unsigned int dbits = (unsigned int)dReg[g * 4L] | ((unsigned int)dReg[g * 4L + 1] << 8); + unsigned int dbits = (unsigned int)(*(const unsigned short*)(dReg + g * 4L)); // #204 review: dReg 16-B aligned, g*4 aligned → single 16-bit d load float d = sharpi_fp16_to_fp32(dbits); float sc0 = d * (float)s[ 0 + isc];