diff --git a/ggml/src/ggml-cuda/common.cuh b/ggml/src/ggml-cuda/common.cuh index ef929d3d7842..4ae9a06bf5ea 100644 --- a/ggml/src/ggml-cuda/common.cuh +++ b/ggml/src/ggml-cuda/common.cuh @@ -173,6 +173,19 @@ static int ggml_cuda_highest_compiled_arch(const int arch) { // --------------------------------------------------------------------------------------------------------- +// GGML_CUDA_BATCH_INVARIANT=1: for batches of 1 to 4 columns, pick kernels whose per-column +// arithmetic does not depend on the column count on the paths this flag covers: the F16 and BF16 +// mat-vec paths it selects, the PTQ1_0 mat-vec (1 to 4 columns; 5 and above take the MMQ tile +// path, which this flag does not touch), and flash attention up to 8 queries. On those paths a +// token decoded alone and a token verified inside a speculative batch see the same logits bit +// for bit. Other weight types and attention shapes outside those kernels can still pick +// batch-dependent kernels, so this is not a whole-model guarantee. Costs some throughput at +// 2 to 4 columns. +static inline bool ggml_cuda_batch_invariant() { + static const bool enabled = getenv("GGML_CUDA_BATCH_INVARIANT") != nullptr; + return enabled; +} + #define MATRIX_ROW_PADDING 512 // last row of quant. matrices is a multiple of this to avoid out-of-bounds memory accesses #define GGML_CUDA_MAX_STREAMS 8 diff --git a/ggml/src/ggml-cuda/fattn-common.cuh b/ggml/src/ggml-cuda/fattn-common.cuh index e67cc7fdf784..ae4318cbc881 100644 --- a/ggml/src/ggml-cuda/fattn-common.cuh +++ b/ggml/src/ggml-cuda/fattn-common.cuh @@ -1154,11 +1154,14 @@ void launch_fattn( // If ntiles_total % blocks_per_wave != 0 then some efficiency is lost due to tail effects. // Test whether parallel_blocks can be set to a higher value for better efficiency. + // Batch-invariant mode: size the KV split as for a single query tile, so the order in which the + // partial softmax results are combined does not depend on how many queries are in the batch. + const int ntiles_dst_eff = ggml_cuda_batch_invariant() ? ntiles_dst / ntiles_x : ntiles_dst; const int blocks_per_wave = nsm * max_blocks_per_sm; int nwaves_best = 0; int efficiency_percent_best = 0; for (int parallel_blocks_test = parallel_blocks; parallel_blocks_test <= ntiles_KV; ++parallel_blocks_test) { - const int nblocks_total = ntiles_dst * parallel_blocks_test; + const int nblocks_total = ntiles_dst_eff * parallel_blocks_test; const int nwaves = (nblocks_total + blocks_per_wave - 1) / blocks_per_wave; const int efficiency_percent = 100 * nblocks_total / (nwaves*blocks_per_wave); diff --git a/ggml/src/ggml-cuda/fattn.cu b/ggml/src/ggml-cuda/fattn.cu index ab7a3b297c07..b5042840da7e 100644 --- a/ggml/src/ggml-cuda/fattn.cu +++ b/ggml/src/ggml-cuda/fattn.cu @@ -460,6 +460,12 @@ static best_fattn_kernel ggml_cuda_get_best_fattn_kernel(const int device, const // If Turing tensor cores are available, use them: if (turing_mma_available(cc) && Q->ne[0] != 40 && Q->ne[0] != 72) { if (can_use_vector_kernel) { + // batch-invariant mode: the same (vector) kernel for up to 8 queries, so a token verified in a + // speculative batch attends with the same arithmetic as a token decoded alone (the whole-model + // guarantee is 1 to 4 queries, see ggml_cuda_batch_invariant) + if (ggml_cuda_batch_invariant() && Q->ne[1] <= 8 && Q->ne[3] == 1) { + return BEST_FATTN_KERNEL_VEC; + } if (!ggml_is_quantized(K->type) && !ggml_is_quantized(V->type)) { if (cc >= GGML_CUDA_CC_ADA_LOVELACE && Q->ne[1] == 1 && Q->ne[3] == 1 && !(gqa_ratio > 4 && K->ne[1] >= 8192)) { return BEST_FATTN_KERNEL_VEC; diff --git a/ggml/src/ggml-cuda/mmvf.cu b/ggml/src/ggml-cuda/mmvf.cu index d7dbc8b99282..b12245e71f63 100644 --- a/ggml/src/ggml-cuda/mmvf.cu +++ b/ggml/src/ggml-cuda/mmvf.cu @@ -821,6 +821,9 @@ bool ggml_cuda_should_use_mmvf(enum ggml_type type, int cc, const int64_t * src0 if (GGML_CUDA_CC_IS_NVIDIA(cc)) { const bool src0_small = (src0_ne[1] <= 512 || src0_ne[2]*src0_ne[3] == 1); if (ampere_mma_available(cc)) { + if (ggml_cuda_batch_invariant()) { + return src0_small && ne11 <= MMVF_MAX_BATCH_SIZE; + } return src0_small && ne11 == 1; } if (cc >= GGML_CUDA_CC_ADA_LOVELACE) { @@ -847,6 +850,11 @@ bool ggml_cuda_should_use_mmvf(enum ggml_type type, int cc, const int64_t * src0 if (GGML_CUDA_CC_IS_NVIDIA(cc)) { const bool src0_small = (src0_ne[1] <= 512 || src0_ne[2]*src0_ne[3] == 1); if (ampere_mma_available(cc)) { + // a few dozen rows (the qwen35 gated-delta-net gate projections) run faster as a + // mat-vec than through the tensor-core path at 2 to 8 columns: 3.4 vs 10.5 us on an RTX 3060 + if (ggml_cuda_batch_invariant() || src0_ne[1] <= 64) { + return src0_small && ne11 <= MMVF_MAX_BATCH_SIZE; + } return src0_small && ne11 == 1; } if (cc >= GGML_CUDA_CC_ADA_LOVELACE) { diff --git a/ggml/src/ggml-cuda/mmvq-ptq1_0.cuh b/ggml/src/ggml-cuda/mmvq-ptq1_0.cuh new file mode 100644 index 000000000000..baa90977b0ef --- /dev/null +++ b/ggml/src/ggml-cuda/mmvq-ptq1_0.cuh @@ -0,0 +1,465 @@ +// PTQ1_0 mat-vec inner loop on a planar-transposed Q8_1 activation layout. +// +// Why: the stock mmvq path hands each thread one 128-weight PTQ1_0 block and +// walks the activations as 32 scattered 4-byte loads per column out of +// 36-byte block_q8_1 structs. The weight decode is shared across columns but +// the activation traffic is not, so every extra column costs a full pass +// (measured 2.2x at 2 columns, 3.3x at 3 on an RTX 3060). Here the activations +// are stored so that the 128 quants a thread needs are 8 aligned 16-byte +// pieces, one per plane, and the 4 (d, s) scales are one more 16-byte piece. +// Adjacent threads read adjacent 16-byte pieces of the same plane, so a warp +// load touches 4 cache lines instead of 32, and a thread reuses each piece +// for every row it owns. The decode, the dp4a sequence and the fp32 epilogue +// are the same operations in the same order as the single-column kernel, so +// every column count produces the same bits for a given column. +// +// PT layout, per activation column (all sizes for the padded row length): +// plane t (t = 0..7): nblk * 16 bytes, byte b of block kb is the quant of +// element kb*128 + t*16 + b +// plane 8: nblk * 16 bytes, block kb holds 4 half2 (d, s), one +// per 32-element sub-block +// Column stride is 9 * nblk * 16 = padded_row * 9/8 bytes, exactly the +// block_q8_1 stride (padded_row/32 blocks of 36 bytes), so every stride the +// mmvq launcher computes in block_q8_1 units stays valid. +#pragma once + +#include "common.cuh" +#include "unary.cuh" +#include "vecdotq.cuh" + +#define PTQ1_0_PT_PLANES 9 + +// dedicated 2D kernel geometry, see mul_mat_vec_ptq1_0_pt below +#define PTQ1_0_PT_THREADS 128 +#define PTQ1_0_PT_MAX_ROWS 16 +#define PTQ1_0_PT_MAX_COLS 4 // 1 to 4 columns here; 5 and above take the MMQ tile path (ggml_cuda_should_use_mmvq) +#define PTQ1_0_PT_SMEM_FLOATS 4096 // 16 KiB of partial sums per weight matrix: the target when choosing rows per CTA + +// the PT path is CUDA only; HIP keeps the block_q8_1 layout and the old vec_dot +static constexpr __host__ __device__ bool ptq1_0_pt_enabled() { +#if defined(GGML_USE_HIP) + return false; +#else + return true; +#endif +} + +// rows of the weight matrix one thread handles per K block: more rows reuse +// each activation piece more often, fewer rows keep the register count down +static constexpr __host__ __device__ int ptq1_0_pt_rows_per_block(const int ncols_dst) { + return ncols_dst == 1 ? 1 : 2; +} + +// number of 128-element blocks in a row padded to MATRIX_ROW_PADDING +static __host__ __device__ __forceinline__ int ptq1_0_pt_nblk(const int ncols_x) { + return ((ncols_x + MATRIX_ROW_PADDING - 1) / MATRIX_ROW_PADDING) * (MATRIX_ROW_PADDING / QK_PTQ1_0); +} + +static __device__ __forceinline__ int int4_at(const int4 & v, const int k) { + switch (k & 3) { + case 0: return v.x; + case 1: return v.y; + case 2: return v.z; + default: return v.w; + } +} + +// four trits (0, 1, 2) packed as bytes -> the weights (-1, 0, 1) as signed bytes. +// t + 127 never carries out of its byte, and flipping the top bit maps +// 127, 128, 129 to -1, 0, 1: two integer ops instead of a byte-wise subtract. +static __device__ __forceinline__ int ptq1_0_trits_to_weights(const int q) { + return (int) (((uint32_t) q + 0x7F7F7F7Fu) ^ 0x80808080u); +} + +// one base-3 digit step on four bytes held as two 16-bit-lane words: +// returns the weights (-1, 0, 1) as four signed bytes, advances the remainders +static __device__ __forceinline__ int ptq1_0_trit_step(uint32_t & vlo, uint32_t & vhi) { + const uint32_t wlo = vlo * 3; + const uint32_t whi = vhi * 3; + vlo = wlo & 0x00FF00FF; + vhi = whi & 0x00FF00FF; + return ptq1_0_trits_to_weights(__byte_perm(wlo, whi, 0x7531)); +} + +// Dot products of nrows PTQ1_0 blocks with the same block index of ncols +// activation columns. bq[i] points at the weight block of row i, ycol[j] at +// the PT column base of column j, kbx is the block index along K. +// +// The integer sum of each 32-element sub-block k is folded into the fp32 +// accumulator as soon as the sub-block is complete, in the order k = 0..3, +// which is the expression acc = sum_k d8_k * sumi_k of the block_q8_1 kernel. +template +static __device__ __forceinline__ void ptq1_0_pt_block_dot( + const block_ptq1_0 * const (&bq)[nrows], + const char * const (&ycol)[ncols], + const int kbx, const int nblk, + float (&result)[ncols][nrows]) { + int4 dsraw[ncols]; +#pragma unroll + for (int j = 0; j < ncols; ++j) { + dsraw[j] = *((const int4 *) ycol[j] + 8*nblk + kbx); + } + + int sumi[ncols][nrows]; + float acc[ncols][nrows]; +#pragma unroll + for (int j = 0; j < ncols; ++j) { +#pragma unroll + for (int i = 0; i < nrows; ++i) { + sumi[j][i] = 0; + acc[j][i] = 0.0f; + } + } + + auto fold = [&](const int k) { +#pragma unroll + for (int j = 0; j < ncols; ++j) { + const float d8 = __low2float(((const half2 *) &dsraw[j])[k]); +#pragma unroll + for (int i = 0; i < nrows; ++i) { + acc[j][i] = __fmaf_rn(d8, (float) sumi[j][i], acc[j][i]); // one FFMA in every instantiation + sumi[j][i] = 0; + } + } + }; + + // qs[0..15]: four groups of four bytes, five trits each: element 16*t + 4*g + b, + // plane t holds words 4*t .. 4*t+3 + uint32_t vlo[nrows][4]; + uint32_t vhi[nrows][4]; +#pragma unroll + for (int i = 0; i < nrows; ++i) { +#pragma unroll + for (int g = 0; g < 4; ++g) { + const uint32_t packed = get_int_b4(bq[i]->qs, g); + vlo[i][g] = __byte_perm(packed, 0, 0x4140); + vhi[i][g] = __byte_perm(packed, 0, 0x4342); + } + } +#pragma unroll + for (int t = 0; t < 5; ++t) { + int4 u[ncols]; +#pragma unroll + for (int j = 0; j < ncols; ++j) { + u[j] = *((const int4 *) ycol[j] + t*nblk + kbx); + } +#pragma unroll + for (int i = 0; i < nrows; ++i) { +#pragma unroll + for (int g = 0; g < 4; ++g) { + const int q = ptq1_0_trit_step(vlo[i][g], vhi[i][g]); +#pragma unroll + for (int j = 0; j < ncols; ++j) { + sumi[j][i] = ggml_cuda_dp4a(q, int4_at(u[j], g), sumi[j][i]); + } + } + } + if (t == 1) { + fold(0); // elements 0..31 done + } + if (t == 3) { + fold(1); // elements 32..63 done + } + } + + // qs[16..23]: two groups of four bytes, five trits each: element 80 + 8*t + 4*g + b, + // words 20..29 live in planes 5, 6 and the lower half of 7 + uint32_t vlo2[nrows][2]; + uint32_t vhi2[nrows][2]; +#pragma unroll + for (int i = 0; i < nrows; ++i) { +#pragma unroll + for (int g = 0; g < 2; ++g) { + const uint32_t packed = get_int_b4(bq[i]->qs + 16, g); + vlo2[i][g] = __byte_perm(packed, 0, 0x4140); + vhi2[i][g] = __byte_perm(packed, 0, 0x4342); + } + } + int4 u2[ncols]; +#pragma unroll + for (int t = 0; t < 5; ++t) { + if (t % 2 == 0) { +#pragma unroll + for (int j = 0; j < ncols; ++j) { + u2[j] = *((const int4 *) ycol[j] + (5 + t/2)*nblk + kbx); + } + } +#pragma unroll + for (int i = 0; i < nrows; ++i) { +#pragma unroll + for (int g = 0; g < 2; ++g) { + const int q = ptq1_0_trit_step(vlo2[i][g], vhi2[i][g]); + const int w = 20 + 2*t + g; // word index within the 128-element block +#pragma unroll + for (int j = 0; j < ncols; ++j) { + sumi[j][i] = ggml_cuda_dp4a(q, int4_at(u2[j], w & 3), sumi[j][i]); + } + } + } + if (t == 1) { + fold(2); // elements 64..95 done (words 16..23) + } + } + + // qh: two bytes, four trits each, interleaved: element 120 + 2*t + h -> words 30, 31, + // the upper half of plane 7 that u2 still holds +#pragma unroll + for (int i = 0; i < nrows; ++i) { + uint32_t v = (uint32_t) bq[i]->qh[0] | ((uint32_t) bq[i]->qh[1] << 16); +#pragma unroll + for (int t = 0; t < 4; t += 2) { + const uint32_t w0 = v * 3; + v = w0 & 0x00FF00FF; + const uint32_t w1 = v * 3; + v = w1 & 0x00FF00FF; + const int q = ptq1_0_trits_to_weights(__byte_perm(w0, w1, 0x7531)); +#pragma unroll + for (int j = 0; j < ncols; ++j) { + sumi[j][i] = ggml_cuda_dp4a(q, int4_at(u2[j], 2 + t/2), sumi[j][i]); + } + } + } + fold(3); // elements 96..127 done + +#pragma unroll + for (int j = 0; j < ncols; ++j) { +#pragma unroll + for (int i = 0; i < nrows; ++i) { + result[j][i] = __fmul_rn((float) bq[i]->d, acc[j][i]); + } + } +} + +// --------------------------------------------------------------------------- +// Dedicated PTQ1_0 mat-vec for plain 2D MUL_MAT (no batch dims, no expert ids). +// +// The generic mmvq kernel gives every thread of a 128-thread block one +// 128-weight K block of the same row, so a K = 5120 projection (40 blocks per +// row) keeps 40 of 128 threads busy and K = 17408 (136 blocks) keeps 53%. Here +// the work items are (row group, K block) pairs of `rows_per_cta` rows +// flattened into one index space, with rows_per_cta chosen on the host so that +// the items fill whole 128-thread iterations where possible. A thread handles +// ROWS adjacent rows per item so that each activation piece it loads serves +// ROWS rows. Each thread writes one fp32 partial per (row, column, K block) to +// shared memory and one warp per (row, column) sums them in a fixed order +// (lane-strided sequential, then a butterfly). That order depends only on the +// weight shape, so the result for a column is the same bits for every column +// count. +// --------------------------------------------------------------------------- + +// rows_per_cta: fill whole 128-thread iterations where possible, within the shared memory budget +static __host__ int ptq1_0_pt_rows_per_cta(const int blocks_per_row, const int ncols_dst, const int nrows_x, const int rows_per_item) { + int rmax = PTQ1_0_PT_SMEM_FLOATS / (ncols_dst * blocks_per_row); + rmax = rmax < rows_per_item ? rows_per_item : (rmax > PTQ1_0_PT_MAX_ROWS ? PTQ1_0_PT_MAX_ROWS : rmax); + rmax -= rmax % rows_per_item; + int best = rows_per_item; + double best_util = 0.0; + for (int r = rows_per_item; r <= rmax; r += rows_per_item) { + const int items = (r / rows_per_item) * blocks_per_row; + const int iters = (items + PTQ1_0_PT_THREADS - 1) / PTQ1_0_PT_THREADS; + const double util = (double) items / (double) (iters * PTQ1_0_PT_THREADS); + if (util > best_util + 1e-9) { + best_util = util; + best = r; + } + if (util > 0.999) { + break; + } + } + GGML_UNUSED(nrows_x); + return best; +} + +// rows per work item: 4 (independent blocks per thread, activation reuse across rows); the 2-row branch served +// 5 to 8 columns before the column limit dropped to 4 and is kept only so the arithmetic stays in one place +static constexpr __host__ __device__ int ptq1_0_pt_rows_per_item(const int ncols_dst) { + return ncols_dst <= 4 ? 4 : 2; +} + +// dynamic shared memory the launch requests for a shape: one fp32 partial per (column, row, K block) for the +// rows of one CTA, twice that with gate fusion. The entry guard and the launcher both call this, so the guard +// budgets exactly what the launch asks for; rows_per_cta is never below rows_per_item, which is what lets a +// long enough K exceed the PTQ1_0_PT_SMEM_FLOATS target. +static __host__ size_t ptq1_0_pt_smem_bytes(const int blocks_per_row, const int ncols_dst, const int nrows_x, const bool has_gate) { + const int rows_per_cta = ptq1_0_pt_rows_per_cta(blocks_per_row, ncols_dst, nrows_x, ptq1_0_pt_rows_per_item(ncols_dst)); + return (size_t) ncols_dst * rows_per_cta * blocks_per_row * sizeof(float) * (has_gate ? 2 : 1); +} + +template +__launch_bounds__(PTQ1_0_PT_THREADS, (ncols <= 2 ? 4 : (ncols <= 4 ? 3 : 2))) +static __global__ void mul_mat_vec_ptq1_0_pt( + const void * vx_ptr, const void * vy_ptr, const ggml_cuda_mm_fusion_args_device fusion, + float * dst_ptr, + const int ncols_x, const int nrows_x, const int stride_row_x, const int stride_col_y, const int stride_col_dst, + const int rows_per_cta, const uint3 bpr_fd) { + // GGML_CUDA_RESTRICT stays off the formal parameters: it expands differently in the host pass and in + // the Hopper-or-newer device pass with PDL, and the generated host stub then fails to match the + // template. Same pattern as mul_mat_vec_q. + const void * GGML_CUDA_RESTRICT vx = vx_ptr; + const void * GGML_CUDA_RESTRICT vy = vy_ptr; + float * GGML_CUDA_RESTRICT dst = dst_ptr; + extern __shared__ float partials_dyn[]; + float * partials = partials_dyn; // [ncols][rows_per_cta][bpr] + [[maybe_unused]] float * partials_gate = partials_dyn + ncols*rows_per_cta*(ncols_x / QK_PTQ1_0); + + const int bpr = ncols_x / QK_PTQ1_0; // K blocks per row + const int nblk = ptq1_0_pt_nblk(ncols_x); // plane stride of the PT layout + const int row0 = rows_per_cta * blockIdx.x; + const int tid = threadIdx.x; + + const char * ycol[ncols]; +#pragma unroll + for (int j = 0; j < ncols; ++j) { + ycol[j] = (const char *) ((const block_q8_1 *) vy + j*stride_col_y); + } + + const int n_items = (rows_per_cta / ROWS) * bpr; + for (int idx = tid; idx < n_items; idx += PTQ1_0_PT_THREADS) { + const int rg = fastdiv((uint32_t) idx, bpr_fd); // row group within the CTA + const int kbx = idx - rg*bpr; + + const block_ptq1_0 * bq[ROWS]; +#pragma unroll + for (int i = 0; i < ROWS; ++i) { + int row = row0 + rg*ROWS + i; + row = row < nrows_x ? row : nrows_x - 1; // clamp the tail, that result is not written + bq[i] = (const block_ptq1_0 *) vx + (int64_t) row*stride_row_x + kbx; + } + float dots[ncols][ROWS]; + ptq1_0_pt_block_dot(bq, ycol, kbx, nblk, dots); +#pragma unroll + for (int j = 0; j < ncols; ++j) { +#pragma unroll + for (int i = 0; i < ROWS; ++i) { + partials[(j*rows_per_cta + rg*ROWS + i)*bpr + kbx] = dots[j][i]; + } + } + if constexpr (has_gate) { + const block_ptq1_0 * bg[ROWS]; +#pragma unroll + for (int i = 0; i < ROWS; ++i) { + int row = row0 + rg*ROWS + i; + row = row < nrows_x ? row : nrows_x - 1; + bg[i] = (const block_ptq1_0 *) fusion.gate + (int64_t) row*stride_row_x + kbx; + } + ptq1_0_pt_block_dot(bg, ycol, kbx, nblk, dots); +#pragma unroll + for (int j = 0; j < ncols; ++j) { +#pragma unroll + for (int i = 0; i < ROWS; ++i) { + partials_gate[(j*rows_per_cta + rg*ROWS + i)*bpr + kbx] = dots[j][i]; + } + } + } + } + + __syncthreads(); + + // one warp per (row, column), lane-strided sequential sum then butterfly: a fixed order + const int warp = tid / WARP_SIZE; + const int lane = tid % WARP_SIZE; + for (int w = warp; w < rows_per_cta*ncols; w += PTQ1_0_PT_THREADS / WARP_SIZE) { + const int j = w / rows_per_cta; + const int r = w - j*rows_per_cta; + const int row = row0 + r; + + float sum = 0.0f; + [[maybe_unused]] float sum_gate = 0.0f; + for (int kbx = lane; kbx < bpr; kbx += WARP_SIZE) { + sum += partials[(j*rows_per_cta + r)*bpr + kbx]; + if constexpr (has_gate) { + sum_gate += partials_gate[(j*rows_per_cta + r)*bpr + kbx]; + } + } + sum = warp_reduce_sum(sum); + if constexpr (has_gate) { + sum_gate = warp_reduce_sum(sum_gate); + } + + if (lane == 0 && row < nrows_x) { + float result = sum; + if constexpr (has_fusion) { + if (fusion.x_bias) { + result += ((const float *) fusion.x_bias)[j*stride_col_dst + row]; + } + if constexpr (has_gate) { + float gate_value = sum_gate; + if (fusion.gate_bias) { + gate_value += ((const float *) fusion.gate_bias)[j*stride_col_dst + row]; + } + switch (fusion.glu_op) { + case GGML_GLU_OP_SWIGLU: + result *= ggml_cuda_op_silu_single(gate_value); + break; + case GGML_GLU_OP_GEGLU: + result *= ggml_cuda_op_gelu_single(gate_value); + break; + case GGML_GLU_OP_SWIGLU_OAI: + result = ggml_cuda_op_swiglu_oai_single(gate_value, result); + break; + default: + result = result * gate_value; + break; + } + } + } + dst[j*stride_col_dst + row] = result; + } + } +} + +template +static void mul_mat_vec_ptq1_0_pt_launch( + const void * vx, const void * vy, const ggml_cuda_mm_fusion_args_device & fusion, float * dst, + const int ncols_x, const int nrows_x, const int stride_row_x, const int stride_col_y, const int stride_col_dst, + cudaStream_t stream) { + constexpr int ROWS = ptq1_0_pt_rows_per_item(ncols); // independent blocks per thread for latency hiding, activation reuse across rows; 8 spills + const int bpr = ncols_x / QK_PTQ1_0; + const int rows_per_cta = ptq1_0_pt_rows_per_cta(bpr, ncols, nrows_x, ROWS); + const uint3 bpr_fd = init_fastdiv_values((uint32_t) bpr); + const dim3 block_nums((nrows_x + rows_per_cta - 1) / rows_per_cta, 1, 1); + const dim3 block_dims(PTQ1_0_PT_THREADS, 1, 1); + + const bool has_fusion = fusion.gate != nullptr || fusion.x_bias != nullptr || fusion.gate_bias != nullptr; + const size_t smem = ptq1_0_pt_smem_bytes(bpr, ncols, nrows_x, fusion.gate != nullptr); + const ggml_cuda_kernel_launch_params lp = ggml_cuda_kernel_launch_params(block_nums, block_dims, smem, stream); + if (has_fusion) { + GGML_ASSERT(ncols == 1 && "fusion only supported for ncols_dst=1"); + if (fusion.gate != nullptr) { + ggml_cuda_kernel_launch(mul_mat_vec_ptq1_0_pt, lp, + vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, rows_per_cta, bpr_fd); + } else { + ggml_cuda_kernel_launch(mul_mat_vec_ptq1_0_pt, lp, + vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, rows_per_cta, bpr_fd); + } + return; + } + ggml_cuda_kernel_launch(mul_mat_vec_ptq1_0_pt, lp, + vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, rows_per_cta, bpr_fd); +} + +// true when the dedicated kernel handles this call (plain 2D, K a multiple of 128, 1 to 4 columns) +static bool mul_mat_vec_ptq1_0_pt_switch( + const void * vx, const void * vy, const ggml_cuda_mm_fusion_args_device & fusion, float * dst, + const int ncols_x, const int nrows_x, const int ncols_dst, + const int stride_row_x, const int stride_col_y, const int stride_col_dst, + const int nchannels_dst, const int nsamples_dst, cudaStream_t stream) { + if (!ptq1_0_pt_enabled() || nchannels_dst != 1 || nsamples_dst != 1 || ncols_x % QK_PTQ1_0 != 0 || + ncols_dst < 1 || ncols_dst > PTQ1_0_PT_MAX_COLS) { + return false; + } + // the kernel never opts in to more than the default dynamic shared memory per block (48 KiB on the + // supported cards), so a shape whose launch would request more goes to the generic kernel + const size_t smem = ptq1_0_pt_smem_bytes(ncols_x / QK_PTQ1_0, ncols_dst, nrows_x, fusion.gate != nullptr); + if (smem > ggml_cuda_info().devices[ggml_cuda_get_device()].smpb) { + return false; + } + switch (ncols_dst) { + case 1: mul_mat_vec_ptq1_0_pt_launch<1>(vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, stream); break; + case 2: mul_mat_vec_ptq1_0_pt_launch<2>(vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, stream); break; + case 3: mul_mat_vec_ptq1_0_pt_launch<3>(vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, stream); break; + case 4: mul_mat_vec_ptq1_0_pt_launch<4>(vx, vy, fusion, dst, ncols_x, nrows_x, stride_row_x, stride_col_y, stride_col_dst, stream); break; + default: return false; + } + return true; +} diff --git a/ggml/src/ggml-cuda/mmvq.cu b/ggml/src/ggml-cuda/mmvq.cu index bf51b61e17b1..3678c788dd5f 100644 --- a/ggml/src/ggml-cuda/mmvq.cu +++ b/ggml/src/ggml-cuda/mmvq.cu @@ -1,4 +1,5 @@ #include "mmvq.cuh" +#include "mmvq-ptq1_0.cuh" #include "quantize.cuh" #include "unary.cuh" #include "vecdotq.cuh" @@ -296,7 +297,11 @@ bool ggml_cuda_should_use_mmvq(enum ggml_type type, int cc, int64_t ne11) { } #if !defined(GGML_USE_HIP) if (type == GGML_TYPE_PTQ1_0 && GGML_CUDA_CC_IS_NVIDIA(cc) && cc >= GGML_CUDA_CC_TURING) { - return ne11 <= 7; + // the PT mat-vec path shares the weight decode across columns and stays ahead of the + // MMQ tile path up to 4 columns; from 5 on the branch-free MMQ tile loader is faster + // (RTX 3060, K = 5120 shapes: a batch of 8 in 64.4 ms through MMQ against 112.7 ms + // through the 8-column mat-vec, llama-bench pp8), so 5 and above take that path + return ne11 <= PTQ1_0_PT_MAX_COLS; } #endif // k-quants cost more to decode and mvq redoes that per column, so MMQ wins sooner. @@ -404,6 +409,12 @@ static constexpr __device__ int get_mmvq_mmid_max_batch_for_device() { } static constexpr __host__ __device__ int calc_nwarps(ggml_type type, int ncols_dst, mmvq_parameter_table_id table_id, bool small_k = false, bool halve_iters = false) { + if (ptq1_0_pt_enabled() && type == GGML_TYPE_PTQ1_0 && + (table_id == MMVQ_PARAMETERS_GENERIC || table_id == MMVQ_PARAMETERS_TURING)) { + // one K partition for every column count, so the fp32 sums of a column + // do not depend on how many columns share the launch (batch invariance) + return ncols_dst <= MMVQ_MAX_BATCH_SIZE ? 4 : 1; + } if (table_id == MMVQ_PARAMETERS_GENERIC) { switch (ncols_dst) { case 1: @@ -534,7 +545,11 @@ static constexpr __host__ __device__ int calc_nwarps(ggml_type type, int ncols_d return 1; } -static constexpr __host__ __device__ int calc_rows_per_block(int ncols_dst, int table_id, bool small_k = false, int nwarps = 1) { +static constexpr __host__ __device__ int calc_rows_per_block(ggml_type type, int ncols_dst, int table_id, bool small_k = false, int nwarps = 1) { + if (ptq1_0_pt_enabled() && type == GGML_TYPE_PTQ1_0 && + (table_id == MMVQ_PARAMETERS_GENERIC || table_id == MMVQ_PARAMETERS_TURING)) { + return ptq1_0_pt_rows_per_block(ncols_dst); + } if (table_id == MMVQ_PARAMETERS_GENERIC || table_id == MMVQ_PARAMETERS_GCN || table_id == MMVQ_PARAMETERS_TURING || table_id == MMVQ_PARAMETERS_GB10) { switch (ncols_dst) { case 1: @@ -573,7 +588,7 @@ static __global__ void mul_mat_vec_q( constexpr int vdr = get_vdr_mmvq(type); constexpr mmvq_parameter_table_id table_id = get_device_table_id(); constexpr int nwarps = calc_nwarps(type, ncols_dst, table_id, small_k, halve_iters); - constexpr int rows_per_cuda_block = calc_rows_per_block(ncols_dst, table_id, small_k, nwarps); + constexpr int rows_per_cuda_block = calc_rows_per_block(type, ncols_dst, table_id, small_k, nwarps); constexpr int warp_size = ggml_cuda_get_physical_warp_size(); constexpr vec_dot_q_cuda_t vec_dot_q_cuda = get_vec_dot_q_cuda(type); @@ -716,28 +731,49 @@ static __global__ void mul_mat_vec_q( const int kqs = vdr * (tid % (qi/vdr)); #if !defined(GGML_USE_HIP) - if constexpr (type == GGML_TYPE_PTQ1_0 && ncols_dst > 1 && ncols_dst <= 3) { + if constexpr (type == GGML_TYPE_PTQ1_0) { + // activations arrive in the PT layout (see mmvq-ptq1_0.cuh): every + // column count runs this same code, one thread per 128-weight block + const int nblk = ptq1_0_pt_nblk(ncols_x); + const char * ycol[ncols_dst]; +# pragma unroll + for (int j = 0; j < ncols_dst; ++j) { + ycol[j] = (const char *) (y + j*stride_col_y); + } + const block_ptq1_0 * bq[rows_per_cuda_block]; # pragma unroll for (int i = 0; i < rows_per_cuda_block; ++i) { - float dots[ncols_dst]; - vec_dot_ptq1_0_q8_1_multi(vx, &y[kby], kbx_offset + i * stride_row_x + kbx, kqs, - stride_col_y, dots); + bq[i] = (const block_ptq1_0 *) vx + kbx_offset + i*stride_row_x + kbx; + } + float dots[ncols_dst][rows_per_cuda_block]; + ptq1_0_pt_block_dot(bq, ycol, kbx, nblk, dots); # pragma unroll - for (int j = 0; j < ncols_dst; ++j) { - tmp[j][i] += dots[j]; + for (int j = 0; j < ncols_dst; ++j) { +# pragma unroll + for (int i = 0; i < rows_per_cuda_block; ++i) { + tmp[j][i] += dots[j][i]; } + } - if constexpr (has_fusion) { - if constexpr (has_gate) { - vec_dot_ptq1_0_q8_1_multi(vgate, &y[kby], kbx_offset + i * stride_row_x + kbx, kqs, - stride_col_y, dots); + if constexpr (has_fusion) { + if constexpr (has_gate) { + const block_ptq1_0 * bg[rows_per_cuda_block]; +# pragma unroll + for (int i = 0; i < rows_per_cuda_block; ++i) { + bg[i] = (const block_ptq1_0 *) vgate + kbx_offset + i*stride_row_x + kbx; + } + ptq1_0_pt_block_dot(bg, ycol, kbx, nblk, dots); # pragma unroll - for (int j = 0; j < ncols_dst; ++j) { - tmp_gate[j][i] += dots[j]; + for (int j = 0; j < ncols_dst; ++j) { +# pragma unroll + for (int i = 0; i < rows_per_cuda_block; ++i) { + tmp_gate[j][i] += dots[j][i]; } } } } + GGML_UNUSED(kqs); + GGML_UNUSED(kby); } else #endif { @@ -897,9 +933,31 @@ static __global__ void mul_mat_vec_q_moe( const int kby = kbx * (qk/QK8_1); const int kqs = vdr * (threadIdx.x % (qi/vdr)); +#if !defined(GGML_USE_HIP) + if constexpr (type == GGML_TYPE_PTQ1_0) { + // PT activation layout, see mmvq-ptq1_0.cuh + const int nblk = ptq1_0_pt_nblk(ncols_x); + const char * ycol[1] = { (const char *) y }; + const block_ptq1_0 * bq[c_rows_per_block]; +#pragma unroll + for (int i = 0; i < c_rows_per_block; ++i) { + bq[i] = (const block_ptq1_0 *) vx + kbx_offset + i*stride_row_x + kbx; + } + float dots[1][c_rows_per_block]; + ptq1_0_pt_block_dot<1, c_rows_per_block>(bq, ycol, kbx, nblk, dots); +#pragma unroll + for (int i = 0; i < c_rows_per_block; ++i) { + tmp[i] += dots[0][i]; + } + GGML_UNUSED(kqs); + GGML_UNUSED(kby); + } else +#endif + { #pragma unroll - for (int i = 0; i < c_rows_per_block; ++i) { - tmp[i] += vec_dot_q_cuda(vx, &y[kby], kbx_offset + i*stride_row_x + kbx, kqs); + for (int i = 0; i < c_rows_per_block; ++i) { + tmp[i] += vec_dot_q_cuda(vx, &y[kby], kbx_offset + i*stride_row_x + kbx, kqs); + } } } @@ -922,7 +980,7 @@ static std::pair calc_launch_params( const int ncols_dst, const int nrows_x, const int nchannels_dst, const int nsamples_or_ntokens, const int warp_size, const mmvq_parameter_table_id table_id, const bool small_k = false, const bool halve_iters = false) { const int nwarps = calc_nwarps(type, ncols_dst, table_id, small_k, halve_iters); - const int rpb = calc_rows_per_block(ncols_dst, table_id, small_k, nwarps); + const int rpb = calc_rows_per_block(type, ncols_dst, table_id, small_k, nwarps); const int64_t nblocks = (nrows_x + rpb - 1) / rpb; const dim3 block_nums(nblocks, nchannels_dst, nsamples_or_ntokens); const dim3 block_dims(warp_size, nwarps, 1); @@ -1003,6 +1061,17 @@ static void mul_mat_vec_q_switch_ncols_dst( GGML_ASSERT(ncols_x % ggml_blck_size(type) == 0); GGML_ASSERT(ncols_dst <= MMVQ_MAX_BATCH_SIZE); +#if !defined(GGML_USE_HIP) + if constexpr (type == GGML_TYPE_PTQ1_0) { + // plain 2D PTQ1_0 mat-vec: dedicated kernel with full lane utilization, see mmvq-ptq1_0.cuh + if (!ids && mul_mat_vec_ptq1_0_pt_switch(vx, vy, fusion, dst, ncols_x, nrows_x, ncols_dst, + stride_row_x, stride_col_y, stride_col_dst, + nchannels_dst, nsamples_dst, stream)) { + return; + } + } +#endif + const uint3 nchannels_y_fd = ids ? init_fastdiv_values(nchannels_y) : make_uint3(0, 0, 0); const uint3 channel_ratio_fd = ids ? make_uint3(0, 0, 0) : init_fastdiv_values(nchannels_dst / nchannels_x); const uint3 sample_ratio_fd = init_fastdiv_values(nsamples_dst / nsamples_x); diff --git a/ggml/src/ggml-cuda/quantize.cu b/ggml/src/ggml-cuda/quantize.cu index fbbc6314ab33..5c661286aa0b 100644 --- a/ggml/src/ggml-cuda/quantize.cu +++ b/ggml/src/ggml-cuda/quantize.cu @@ -1,4 +1,5 @@ #include "quantize.cuh" +#include "mmvq-ptq1_0.cuh" #include "unary.cuh" #include @@ -51,6 +52,9 @@ static __device__ __forceinline__ float nvfp4_native_scale_error( #endif // CUDART_VERSION >= 12080 #endif // defined(BLACKWELL_MMA_AVAILABLE) +// pt: write the planar-transposed layout consumed by the PTQ1_0 mat-vec path +// (see mmvq-ptq1_0.cuh) instead of block_q8_1; same quantization, same bytes per row +template __launch_bounds__(CUDA_QUANTIZE_BLOCK_SIZE, 1) static __global__ void quantize_q8_1( const float * x_ptr, void * vy_ptr, @@ -92,6 +96,23 @@ static __global__ void quantize_q8_1( const float d = amax / 127.0f; const int8_t q = amax == 0.0f ? 0 : roundf(xi / d); + if constexpr (pt) { + const int64_t row_cont = (i3*ne2.z + i2) * ne1 + i1; + char * ycol = (char *) vy + row_cont * (ne0 * 9 / 8); // same row stride as block_q8_1 + const int64_t nblk = ne0 / QK_PTQ1_0; + const int64_t kb = i0 / QK_PTQ1_0; + const int e = i0 % QK_PTQ1_0; + ycol[((e / 16)*nblk + kb) * 16 + (e % 16)] = q; + + if (iqs > 0) { + return; + } + + half2 * ds = (half2 *) (ycol + 8*nblk*16) + kb*4 + e / QK8_1; + *ds = make_half2(d, sum); + return; + } + y[ib].qs[iqs] = q; if (iqs > 0) { @@ -648,8 +669,12 @@ void quantize_row_q8_1_cuda( const dim3 num_blocks(block_num_x, ne1, ne2*ne3); const dim3 block_size(CUDA_QUANTIZE_BLOCK_SIZE, 1, 1); const ggml_cuda_kernel_launch_params launch_params = ggml_cuda_kernel_launch_params(num_blocks, block_size, 0, stream); - ggml_cuda_kernel_launch(quantize_q8_1, launch_params, x, vy, ne00, s01, s02, s03, ne0, ne1, ne2_fastdiv); - GGML_UNUSED(type_src0); + if (ptq1_0_pt_enabled() && type_src0 == GGML_TYPE_PTQ1_0) { + GGML_ASSERT(ne0 % QK_PTQ1_0 == 0); + ggml_cuda_kernel_launch(quantize_q8_1, launch_params, x, vy, ne00, s01, s02, s03, ne0, ne1, ne2_fastdiv); + return; + } + ggml_cuda_kernel_launch(quantize_q8_1, launch_params, x, vy, ne00, s01, s02, s03, ne0, ne1, ne2_fastdiv); } void quantize_mmq_q8_1_cuda( diff --git a/tests/test-backend-ops.cpp b/tests/test-backend-ops.cpp index 14f654f966ca..ddde50a7ecdc 100644 --- a/tests/test-backend-ops.cpp +++ b/tests/test-backend-ops.cpp @@ -10174,6 +10174,17 @@ static std::vector> make_test_cases_eval() { } } + // PTQ1_0 mat-vec kernel, shared-memory boundary: the launch asks for ncols * rows_per_cta * (K / 128) * 4 bytes, + // twice that with gate fusion, against the 48 KiB default. With one column and 4 rows per CTA the last K that + // fits is 393216 without a gate and 196608 with one; the next K block over each must take the generic kernel. + for (int64_t k : {393216, 393344}) { + test_cases.emplace_back(new test_mul_mat(GGML_TYPE_PTQ1_0, GGML_TYPE_F32, 64, 1, k, {1, 1}, {1, 1})); + } + for (int64_t k : {196608, 196736}) { + test_cases.emplace_back(new test_mul_mat_vec_fusion(GGML_TYPE_PTQ1_0, GGML_GLU_OP_SWIGLU, 1, 64, k, + false, 1, 1, false, false, true, false, {1, 1})); + } + for (auto gate : {GATING_FUNC_SOFTMAX, GATING_FUNC_SIGMOID, GATING_FUNC_SOFTMAX_WEIGHT, GATING_FUNC_SQRT_SOFTPLUS}) { for (bool with_norm : {false, true}) { for (bool bias_probs : {false, true}) { @@ -10476,6 +10487,23 @@ static std::vector> make_test_cases_perf() { } } + // Ternary Bonsai 2 27B (qwen35): the bf16 gated-delta-net gate projections + for (int bs : {1, 2, 3, 4, 8}) { + test_cases.emplace_back(new test_mul_mat(GGML_TYPE_BF16, GGML_TYPE_F32, 48, bs, 5120, {1, 1}, {1, 1})); // ssm_alpha, ssm_beta + } + + // Ternary Bonsai 2 27B (qwen35, PTQ1_0) projections at speculative-decoding batch sizes + for (int bs : {1, 2, 3, 4, 8}) { + for (ggml_type type_a : {GGML_TYPE_PTQ1_0, GGML_TYPE_Q4_0}) { + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 10240, bs, 5120, {1, 1}, {1, 1})); // attn_qkv + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 6144, bs, 5120, {1, 1}, {1, 1})); // attn_gate + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 5120, bs, 6144, {1, 1}, {1, 1})); // ssm_out + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 17408, bs, 5120, {1, 1}, {1, 1})); // ffn_up, ffn_gate + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 5120, bs, 17408, {1, 1}, {1, 1})); // ffn_down + test_cases.emplace_back(new test_mul_mat(type_a, GGML_TYPE_F32, 12288, bs, 5120, {1, 1}, {1, 1})); // attn_q + } + } + // qwen3-30b-a3b for (int bs : {1, 4, 8, 32, 64, 128, 256, 512}) { for (ggml_type type_a : {GGML_TYPE_F32, GGML_TYPE_F16, GGML_TYPE_Q4_0, GGML_TYPE_Q8_0, GGML_TYPE_Q4_K, GGML_TYPE_Q6_K, GGML_TYPE_IQ2_XS}) {