From 98ea4108097bc345baf3970178f7927a5b77011d Mon Sep 17 00:00:00 2001 From: sudoingX <200180104+sudoingX@users.noreply.github.com> Date: Sat, 19 Sep 2026 07:08:01 +0000 Subject: [PATCH 1/5] Add: planar-transposed activation layout for the PTQ1_0 mat-vec path The mmvq path gave every thread one 128-weight PTQ1_0 block and read the activations as 32 scattered 4-byte loads per column out of 36-byte block_q8_1 structs, so every extra column cost a full pass: 47 / 104 / 154 / 201 us for 1 / 2 / 3 / 4 columns on a 4096 x 14336 matrix (RTX 3060, test-backend-ops perf). Quantize PTQ1_0 activations into a planar-transposed layout instead: per column, 8 planes of 16-byte quant pieces and one plane of the 4 (d, s) half2 scales, 16 bytes per 128-element block each, with the same bytes per row and the same column stride as block_q8_1. A thread then reads its block's activations as 9 aligned 16-byte loads with adjacent lanes on adjacent pieces: 46 / 55 / 65 / 75 us. The trit decode is unchanged; the byte-wise minus one is (q + 0x7F7F7F7F) ^ 0x80808080 instead of __vsub4. nwarps is pinned to 4 for every column count so the K partition, and with it the fp32 summation order, does not depend on the batch. mmvq now handles PTQ1_0 up to 8 columns (MMQ took 395 us at 8). HIP keeps the old layout. --- ggml/src/ggml-cuda/mmvq-ptq1_0.cuh | 224 +++++++++++++++++++++++++++++ ggml/src/ggml-cuda/mmvq.cu | 92 +++++++++--- ggml/src/ggml-cuda/quantize.cu | 29 +++- 3 files changed, 325 insertions(+), 20 deletions(-) create mode 100644 ggml/src/ggml-cuda/mmvq-ptq1_0.cuh 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 000000000..dfe9e2882 --- /dev/null +++ b/ggml/src/ggml-cuda/mmvq-ptq1_0.cuh @@ -0,0 +1,224 @@ +// 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 "vecdotq.cuh" + +#define PTQ1_0_PT_PLANES 9 + +// 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] += d8 * (float) sumi[j][i]; + 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] = (float) bq[i]->d * acc[j][i]; + } + } +} diff --git a/ggml/src/ggml-cuda/mmvq.cu b/ggml/src/ggml-cuda/mmvq.cu index bf51b61e1..5044cd77c 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,9 @@ 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 the full mmvq batch + return ne11 <= MMVQ_MAX_BATCH_SIZE; } #endif // k-quants cost more to decode and mvq redoes that per column, so MMQ wins sooner. @@ -404,6 +407,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 +543,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 +586,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 +729,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) { # pragma unroll - for (int j = 0; j < ncols_dst; ++j) { - tmp_gate[j][i] += dots[j]; + 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 +931,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] += 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] += 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); + } } } @@ -922,7 +978,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); diff --git a/ggml/src/ggml-cuda/quantize.cu b/ggml/src/ggml-cuda/quantize.cu index fbbc6314a..5c661286a 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( -- 2.34.1