From cff66d2065b576d9ce79e440499b0d040acf93e0 Mon Sep 17 00:00:00 2001 From: Lyes SBAHI <143282344+lyesrock@users.noreply.github.com> Date: Sun, 25 Oct 2026 14:55:49 +0000 Subject: [PATCH] DeepSeek engine opt-in: the busy experts of a prompt take FP16 GEMMs, not the tiles (DS4_FAT_ROWS=N, default 0: byte-identical path, the packed tiles run exactly as before). An expert of N rows or more has its weights decoded to FP16 once (iq::f16_rows, the values the tiles multiply), the rows gathered (gather_f16), cublasGemmEx, clamped Swiglu in FP16 (swiglu_f16); SetArgs skip_ge leaves those experts out of the tiles. The tiles decode an expert's weights again for every 64 of its rows; a popular expert of a long chunk has hundreds. Measured on the RTX 5090 with the IST-DAS 2.25 bpw pack, 0.2.7.1 base, DS4_FAT_ROWS=64, ABBA (four fresh boots alternating off/on, a fresh random 13K-token prompt every measurement, slots and page cache cleared, nice=0): the cold 13K prompt 1 144 / 1 145 off against 1 170 / 1 206 on (+3.8%), repeat prompts 1 561-1 594 off against 1 642-1 689 on (+5.7%), decode unchanged (28-33 short, 47-51 warm). The same products, the sums' order differs: greedy output can diverge after the first rounding (wording equivalent); on 12 greedy prompts at temperature 0 the off side is indistinguishable from the stock engine, the differences inside the stock's own run-to-run band. --- include/strata/dsv41/exl3.hpp | 1 + src/dsv41/ds4_gguf.cu | 1 + src/dsv41/ds4_iq.cu | 64 +++++++++++++++++++++ src/dsv41/ds4_iq.cuh | 9 +++ src/dsv41/ds4_model.cu | 101 ++++++++++++++++++++++++++++++++++ 5 files changed, 176 insertions(+) diff --git a/include/strata/dsv41/exl3.hpp b/include/strata/dsv41/exl3.hpp index f11176de193..afd6def60cc 100644 --- a/include/strata/dsv41/exl3.hpp +++ b/include/strata/dsv41/exl3.hpp @@ -158,6 +158,7 @@ struct SetArgs { // their shared-memory tile anyway, so the same values come out; this only saves the F32 read of the same rows by // every column block of an expert (18-36 reads of 1.8-3.7 MB on a prompt) - MarLin's "small operand cheap". const uint16_t* x16 = nullptr; + int skip_ge = 0; // (ds4's low-bit GGUF experts) > 0: experts of at least this many rows are the caller's (FP16 GEMMs) #if defined(STRATA_METAL) // (only on Metal: keeps the CUDA / HIP struct as it was) const int* hbounds = nullptr; // the same bounds on the host: [n_exp + 1], absolute rows int total_rows = 0; // the set's rows in all diff --git a/src/dsv41/ds4_gguf.cu b/src/dsv41/ds4_gguf.cu index 8e1d605be7f..f537a3ce02a 100644 --- a/src/dsv41/ds4_gguf.cu +++ b/src/dsv41/ds4_gguf.cu @@ -1563,6 +1563,7 @@ iq::SetArgs iq_set(const x3::SetArgs& sa) { // a set as ds4_iq's arguments (th q.n_ff = sa.n_ff; q.limit = sa.limit; q.xf16 = sa.x16; + q.skip_ge = sa.skip_ge; #if defined(STRATA_METAL) q.hbounds = sa.hbounds; q.total_rows = sa.total_rows; diff --git a/src/dsv41/ds4_iq.cu b/src/dsv41/ds4_iq.cu index 186aedb9d96..be5d8683701 100644 --- a/src/dsv41/ds4_iq.cu +++ b/src/dsv41/ds4_iq.cu @@ -742,6 +742,7 @@ __global__ void __launch_bounds__(kThreads) set_kernel(const SetArgs a) { const int lo = a.bounds[e] - a.row0, hi = a.bounds[e + 1] - a.row0; const int r0 = lo + blockIdx.x * RM; if (r0 >= hi || hi - lo <= a.skip_le) return; // (uniform; skip_le: the light kernels' experts) + if (a.skip_ge > 0 && hi - lo >= a.skip_ge) return; // (the caller's FP16 GEMMs') load_tab(tab); const int nr = min(RM, hi - r0); const int c0 = blockIdx.y * BN; @@ -833,6 +834,7 @@ __global__ void __launch_bounds__(kThreads) set_tc_kernel(const SetArgs a) { const int lo = a.bounds[e] - a.row0, hi = a.bounds[e + 1] - a.row0; const int r0 = lo + blockIdx.x * BM; if (r0 >= hi || hi - lo <= a.skip_le) return; // (uniform; skip_le: the light kernels' experts) + if (a.skip_ge > 0 && hi - lo >= a.skip_ge) return; // (the caller's FP16 GEMMs') load_tab(tab); const int nr = min(BM, hi - r0); const int c0 = blockIdx.y * BN; @@ -1618,6 +1620,68 @@ bool set_down(const SetArgs& a, cudaStream_t s) { return true; } +// ---------------------------------------------------------------- the busy experts' FP16 weights (f16_rows) +// 8 values a thread: row r's sub-block kb, piece p - deq8's values rounded as the tensor-core tiles round them +template +__global__ void __launch_bounds__(kThreads) f16_rows_kernel(const uint8_t* __restrict__ w, int k, int n, + uint16_t* __restrict__ dst, int ldd) { + __shared__ uint32_t tab[Ty::TAB > 0 ? Ty::TAB : 1]; + load_tab(tab); + __syncthreads(); + const int nq = k / 8; + const size_t rb = row_bytes_t(k); + for (long long i = (long long) blockIdx.x * blockDim.x + threadIdx.x; i < (long long) n * nq; + i += (long long) gridDim.x * blockDim.x) { + const int r = (int) (i / nq), q = (int) (i % nq); + float w8[8]; + deq8(w + (size_t) r * rb, q >> 2, q & 3, tab, w8); + union { + uint4 v; + uint16_t h[8]; + } u; +#pragma unroll + for (int j = 0; j < 8; ++j) u.h[j] = __half_as_ushort(__float2half_rn(w8[j])); + *reinterpret_cast(dst + (size_t) r * ldd + (size_t) q * 8) = u.v; + } +} +__global__ void __launch_bounds__(kThreads) gather_f16_kernel(const uint16_t* __restrict__ x, int ldx, + const int* __restrict__ idx, int n, int k, + uint16_t* __restrict__ dst) { + const int nq = k / 8; + for (long long i = (long long) blockIdx.x * blockDim.x + threadIdx.x; i < (long long) n * nq; + i += (long long) gridDim.x * blockDim.x) { + const int r = (int) (i / nq), q = (int) (i % nq); + *reinterpret_cast(dst + (size_t) r * k + (size_t) q * 8) = + *reinterpret_cast(x + (size_t) idx[r] * ldx + (size_t) q * 8); + } +} +__global__ void __launch_bounds__(kThreads) swiglu_f16_kernel(const float* __restrict__ gu, int ld, int n, int ff, + float limit, uint16_t* __restrict__ h, int ldh) { + for (long long i = (long long) blockIdx.x * blockDim.x + threadIdx.x; i < (long long) n * ff; + i += (long long) gridDim.x * blockDim.x) { + const int r = (int) (i / ff), c = (int) (i % ff); + const float* g = gu + (size_t) r * ld; + h[(size_t) r * ldh + c] = __half_as_ushort(__float2half_rn(swiglu(g[c], g[ff + c], limit))); + } +} +static unsigned grid_for(long long items) { + return (unsigned) std::max(1, std::min(4096, (items + kThreads - 1) / kThreads)); +} +bool f16_rows(int type, const uint8_t* w, int k, int n, uint16_t* dst, int ldd, cudaStream_t s) { + if (k % 256 != 0 || n < 1 || ldd % 8 != 0 || ((uintptr_t) dst & 15) != 0) return false; + const unsigned blocks = grid_for((long long) n * (k / 8)); +#define FR_(T) f16_rows_kernel<<>>(w, k, n, dst, ldd) + DS4_IQ_SWITCH(type, FR_) +#undef FR_ + return cudaGetLastError() == cudaSuccess; +} +void gather_f16(const uint16_t* x, int ldx, const int* idx, int n, int k, uint16_t* dst, cudaStream_t s) { + gather_f16_kernel<<>>(x, ldx, idx, n, k, dst); +} +void swiglu_f16(const float* gu, int ld, int n, int ff, float limit, uint16_t* h, int ldh, cudaStream_t s) { + swiglu_f16_kernel<<>>(gu, ld, n, ff, limit, h, ldh); +} + bool mv(const MvJob* jobs, int n, cudaStream_t s) { #if defined(STRATA_METAL) if (n > 0 && n <= 64 && ::strata::dsv41::metal_host::iq_mv((const void*) jobs, n)) return true; diff --git a/src/dsv41/ds4_iq.cuh b/src/dsv41/ds4_iq.cuh index 6a22576abe4..8cf48e2c1a7 100644 --- a/src/dsv41/ds4_iq.cuh +++ b/src/dsv41/ds4_iq.cuh @@ -95,6 +95,7 @@ struct SetArgs { int ldo = 0; int tg = 0, td = 0; int skip_le = 0; // (the tiles) the experts of at most this many rows are left out: the light kernels took them + int skip_ge = 0; // (the tiles) > 0: the experts of at least this many rows are left out - the caller's GEMMs #if defined(STRATA_METAL) // (only on Metal: the CUDA and HIP kernels take this struct by value, so its size is part of their ABI) const int* hbounds = nullptr; // the bounds on the host: [n_exp + 1], absolute rows int total_rows = 0; // the set's rows in all @@ -106,6 +107,14 @@ bool set_down(const SetArgs& a, cudaStream_t s); /// count of matrix-core launches so far (both do nothing / return -1 where there are no such tiles; hipk_force returns whether this device has them) bool hipk_force(int mode); long long hipk_launches(); +/// The busy experts' path (a prompt's experts of many rows): one expert's matrix - n rows of k values, row_bytes(k) +/// apart - as FP16 rows (ldd apart), the values the tensor-core tiles multiply (their loader's, rounded to FP16 the same +/// way), for the caller's GEMM: the tiles decode an expert's weights again for every 64 of its rows. +bool f16_rows(int type, const uint8_t* w, int k, int n, uint16_t* dst, int ldd, cudaStream_t s); +/// dst[r] = x[idx[r]] for r < n (FP16 rows of k values; x ldx apart, dst k apart) +void gather_f16(const uint16_t* x, int ldx, const int* idx, int n, int k, uint16_t* dst, cudaStream_t s); +/// h[r][c] = FP16 of the tiles' clamped SwiGLU of gu[r][c] (gate) and gu[r][ff + c] (up), F32 rows ld apart +void swiglu_f16(const float* gu, int ld, int n, int ff, float limit, uint16_t* h, int ldh, cudaStream_t s); /// The light form alone (the set's experts of at most 16 rows; the others untouched), for any of these types and MXFP4 /// (tg / td 39: ds4_gguf.cu's experts): false (nothing launched) when it is switched off - the caller's tiles take every /// expert then; true: they take those of more rows (skip_le 16). diff --git a/src/dsv41/ds4_model.cu b/src/dsv41/ds4_model.cu index 7ffa7cc4f9d..8f072397cc8 100644 --- a/src/dsv41/ds4_model.cu +++ b/src/dsv41/ds4_model.cu @@ -612,6 +612,12 @@ struct Chunk { size_t stage_stride = 0; std::unique_ptr readers; cublasHandle_t blas = nullptr; + // the busy experts' path (run_set, DS4_FAT_ROWS): one expert's weights as FP16 (gate and up, then down), the rows + // gathered for its GEMM, the GEMM's gate | up output - fat_r rows a piece; allocated at the first use + uint16_t *fat_w = nullptr, *fat_x = nullptr; + float* fat_gu = nullptr; + int fat_r = 0; + bool fat_off = false; // host-mapped: a chunk's hand-off to the next group (h | pre | topk | cand), two in flight (slot = chunk % 2), each // with its flag, and the next group's members' acknowledgements [slot][member] (32 apart) uint8_t* hand[2] = {nullptr, nullptr}; @@ -3379,6 +3385,78 @@ struct Model::Impl { } return true; } + // ---- the busy experts (run_set): an expert of many rows as FP16 GEMMs. The tiles decode its weights again for + // every 64 of its rows (a prompt's chunk gives a popular expert hundreds); decoded once, to the FP16 values the tiles + // multiply (iq::f16_rows), the library's GEMM takes them at its own speed. DS4_FAT_ROWS=N (default 0, off; N: + // the experts of N rows or more) takes the experts of N rows or more; the same products, only the sums' order differ + // (3x V100, Maya-XS: 64 rows +6-8% on 8K / 32K prompts; KL vs the reference 0.276, the stock path 0.285 - the same; + // RTX 5090, IST-DAS 2.25 bpw: a 13K-token prompt 1,400 -> 1,625 tok/s at 64 rows, decode unchanged) + static int fat_rows() { + static const int n = std::getenv("DS4_FAT_ROWS") ? std::max(0, std::atoi(std::getenv("DS4_FAT_ROWS"))) : 0; + return n; + } + bool fat_ready(Dev& d, Chunk& c, int D, int FF) { +#if defined(STRATA_USE_HIP) + (void) d, (void) c, (void) D, (void) FF; + return false; +#else + if (fat_rows() <= 0 || c.fat_off) return false; + if (c.fat_w != nullptr) return true; + static const int R = std::getenv("DS4_FAT_PIECE") ? std::max(64, std::atoi(std::getenv("DS4_FAT_PIECE"))) : 1024; + if (cudaSetDevice(d.id) != cudaSuccess || cudaMalloc(&c.fat_w, (size_t) 2 * FF * D * 2) != cudaSuccess || + cudaMalloc(&c.fat_x, (size_t) R * D * 2) != cudaSuccess || + cudaMalloc(&c.fat_gu, (size_t) R * 2 * FF * 4) != cudaSuccess) { + cudaGetLastError(); + std::fprintf(stderr, "ds4: CUDA%d: no VRAM for the busy experts' path (%.0f MB): the tiles take them\n", d.id, + ((double) 2 * FF * D * 2 + (double) R * D * 2 + (double) R * 2 * FF * 4) / 1e6); + c.fat_off = true; + return false; + } + c.fat_r = R; + std::fprintf(stderr, "ds4: CUDA%d busy experts (%d rows or more) as FP16 GEMMs, %d rows a piece (%.0f MB)\n", d.id, + fat_rows(), R, ((double) 2 * FF * D * 2 + (double) R * D * 2 + (double) R * 2 * FF * 4) / 1e6); + return true; +#endif + } + // h16[lr0 + r] = FP16 clamped SwiGLU of the expert's gate and up over x16[row_tok[row + r]], r < n + bool fat_gate_up(Chunk& c, const x3::SetArgs& sa, const uint8_t* blob, int row, int n, int a, int D, int FF, + cudaStream_t s, std::string& err) { + const int tg = gg::expert_type(sa.lay, 0); + if (!gg::iq::f16_rows(tg, blob + sa.lay.trellis[0], D, FF, c.fat_w, D, s) || + !gg::iq::f16_rows(tg, blob + sa.lay.trellis[1], D, FF, c.fat_w + (size_t) FF * D, D, s)) { + err = "the busy experts' path: gate / up not decoded"; + return false; + } + const float one = 1.0f, zero = 0.0f; + for (int p0 = 0; p0 < n; p0 += c.fat_r) { + const int np = std::min(c.fat_r, n - p0); + gg::iq::gather_f16(c.x16, D, c.row_tok + row + p0, np, D, c.fat_x, s); + if (cublasGemmEx(c.blas, CUBLAS_OP_T, CUBLAS_OP_N, 2 * FF, np, D, &one, c.fat_w, CUDA_R_16F, D, c.fat_x, + CUDA_R_16F, D, &zero, c.fat_gu, CUDA_R_32F, 2 * FF, CUBLAS_COMPUTE_32F, + CUBLAS_GEMM_DEFAULT_TENSOR_OP) != CUBLAS_STATUS_SUCCESS) { + err = "the busy experts' path: cublasGemmEx (gate / up) failed"; + return false; + } + gg::iq::swiglu_f16(c.fat_gu, 2 * FF, np, FF, sa.limit, c.h16 + (size_t) (row - a + p0) * FF, FF, s); + } + return true; + } + // y[lr0 + r] = down(h16[lr0 + r]), r < n + bool fat_down(Chunk& c, const x3::SetArgs& sa, const uint8_t* blob, int lr0, int n, int D, int FF, cudaStream_t s, + std::string& err) { + if (!gg::iq::f16_rows(gg::expert_type(sa.lay, 2), blob + sa.lay.trellis[2], FF, D, c.fat_w, FF, s)) { + err = "the busy experts' path: down not decoded"; + return false; + } + const float one = 1.0f, zero = 0.0f; + if (cublasGemmEx(c.blas, CUBLAS_OP_T, CUBLAS_OP_N, D, n, FF, &one, c.fat_w, CUDA_R_16F, FF, + c.h16 + (size_t) lr0 * FF, CUDA_R_16F, FF, &zero, c.y + (size_t) lr0 * D, CUDA_R_32F, D, + CUBLAS_COMPUTE_32F, CUBLAS_GEMM_DEFAULT_TENSOR_OP) != CUBLAS_STATUS_SUCCESS) { + err = "the busy experts' path: cublasGemmEx (down) failed"; + return false; + } + return true; + } // y[r][c] = alpha x16[r] . W[c] for F16 weights W [n][k] (the indexer's weights_proj) bool gemm_f16w(Chunk& c, const uint16_t* W, int n, int k, const uint16_t* x16, int ldx, int rows, float alpha, float* y, int ldy, std::string& err) { @@ -4045,6 +4123,21 @@ struct Model::Impl { // gate / up read the layer's F16 copy of the activations (the same rounding the tiles make for their // shared tile): with 18-36 column blocks an expert the F32 rows were read that many times sa.x16 = strata::ds4::gate_up_x16_default(defaults_026()) ? c.x16 : nullptr; + // the busy experts (DS4_FAT_ROWS=N: those of N rows or more): their weights decoded to FP16 once and + // multiplied by the library's GEMM - the tiles decode an expert's weights again for every 64 of its rows + std::vector fat; + if (sa.lay.fmt == gg::kFmtMXFP4 && !gg::mx_experts(sa.lay) && fat_ready(d, c, D, FF)) { + int rest = 0; + for (int i = i0; i < j; ++i) { + const int r = b[i + 1] - b[i]; + if (r >= fat_rows()) fat.push_back(i); + else rest = std::max(rest, r); + } + if (!fat.empty()) { + sa.skip_ge = fat_rows(); + sa.max_rows = rest; + } + } EkMark ek{}; const bool ekp = prof_chunk; if (ekp) { @@ -4063,8 +4156,16 @@ struct Model::Impl { } if (sa.lay.fmt == gg::kFmtMXFP4) { gg::set_gate_up(sa, c.x, D, c.row_tok + a, c.h16, FF, d.st); + for (const int i : fat) + if (!fat_gate_up(c, sa, sa.wbase + (size_t) (i - i0) * stride, b[i], b[i + 1] - b[i], a, D, FF, + d.st, err)) + return false; if (ekp) cudaEventRecord(ek.e1, d.st); gg::set_down(sa, c.h16, FF, c.y, D, d.st); + for (const int i : fat) + if (!fat_down(c, sa, sa.wbase + (size_t) (i - i0) * stride, b[i] - a, b[i + 1] - b[i], D, FF, + d.st, err)) + return false; } else { x3::set_gate_up(sa, c.x, D, c.row_tok + a, c.h16, FF, d.st); if (ekp) cudaEventRecord(ek.e1, d.st);