← Back to C-Kernel-Engine Docs Doxygen Source Documentation
 
Loading...
Searching...
No Matches
gemm_kernels_q6k_q8k.c
Go to the documentation of this file.
1/**
2 * @file gemm_kernels_q6k_q8k.c
3 * @brief Q6_K (weights) x Q8_K (activations) kernels for inference
4 *
5 * CK-ENGINE KERNEL RULES:
6 * =======================
7 * 1. NO malloc/free - memory via bump allocator, pointers passed in
8 * 2. NO OpenMP - parallelization at orchestrator/codegen layer
9 * 3. API must define: inputs, outputs, workspace, and memory layouts
10 * 4. Pure computation - deterministic, no side effects
11 *
12 * After changes: make test && make llamacpp-parity-full
13 *
14 * Implements decode-style matvec/matmul where weights are Q6_K and the
15 * activations are quantized on-the-fly to Q8_K. This is inference-only;
16 * no backward pass is provided here.
17 *
18 * Q6_K Format (256 weights per block):
19 * - d: FP16 super-block scale
20 * - ql: 128 bytes (low 4 bits of each weight)
21 * - qh: 64 bytes (high 2 bits of each weight)
22 * - scales: 16 int8 sub-block scales
23 *
24 * Q8_K Format (256 weights per block):
25 * - d: FP32 scale
26 * - qs: 256 int8 values
27 * - bsums: 16 int16 block sums
28 */
29
30#include <assert.h>
31#include <math.h>
32#include <string.h>
33#include <stdint.h>
34#include <stddef.h>
35#include <stdlib.h>
36
37#include "ckernel_engine.h"
38#include "ckernel_quant.h"
39
40/* Include SIMD headers based on available extensions */
41#if defined(__AVX512F__) || defined(__AVX2__) || defined(__AVX__) || defined(__SSE4_1__) || defined(__SSSE3__)
42#include <immintrin.h>
43#endif
44#if defined(__ARM_NEON) || defined(__aarch64__)
45#include <arm_neon.h>
46#endif
47
48typedef struct {
49 ggml_half d;
50 int8_t scales[16];
51 uint8_t qs[QK_K];
52} block_q6_K_prepared;
53
54_Static_assert(sizeof(block_q6_K_prepared) == 274,
55 "Q6_K prepared-size contract changed");
56
57size_t ck_q6_k_prepared_block_size(void)
58{
59 return sizeof(block_q6_K_prepared);
60}
61
62const char *ck_q6_k_prepared_provider_name(void)
63{
64#if defined(__AVX512F__) && defined(__AVX512BW__) && \
65 defined(__AVX512VNNI__)
66 return "q6_k_prepared_avx512_vnni_exact";
67#elif defined(__AVX2__)
68 return "q6_k_prepared_avx2_exact";
69#else
70 return "q6_k_prepared_unavailable";
71#endif
72}
73
74void ck_q6_k_prepare_weight(const void *src, void *dst, int N, int K)
75{
76 if (!src || !dst || N <= 0 || K <= 0 || (K % QK_K) != 0) return;
77 const block_q6_K *input = (const block_q6_K *)src;
78 block_q6_K_prepared *output = (block_q6_K_prepared *)dst;
79 const size_t blocks = (size_t)N * (size_t)(K / QK_K);
80
81 for (size_t b = 0; b < blocks; ++b) {
82 output[b].d = input[b].d;
83 memcpy(output[b].scales, input[b].scales, sizeof(output[b].scales));
84 for (int n = 0; n < QK_K; n += 128) {
85 const uint8_t *ql = input[b].ql + n / 2;
86 const uint8_t *qh = input[b].qh + n / 4;
87 for (int l = 0; l < 32; ++l) {
88 output[b].qs[n + l + 0] =
89 (uint8_t)((ql[l] & 0x0f) | (((qh[l] >> 0) & 3) << 4));
90 output[b].qs[n + l + 32] =
91 (uint8_t)((ql[l + 32] & 0x0f) | (((qh[l] >> 2) & 3) << 4));
92 output[b].qs[n + l + 64] =
93 (uint8_t)((ql[l] >> 4) | (((qh[l] >> 4) & 3) << 4));
94 output[b].qs[n + l + 96] =
95 (uint8_t)((ql[l + 32] >> 4) | (((qh[l] >> 6) & 3) << 4));
96 }
97 }
98 }
99}
100
101/* Forward declarations for SIMD implementations */
102void gemv_q6_k_q8_k_avx512(float *y, const void *W, const void *x_q8, int M, int K);
103void gemv_q6_k_q8_k_avx512_vbmi(float *y, const void *W, const void *x_q8, int M, int K);
104void gemv_q6_k_q8_k_avx2(float *y, const void *W, const void *x_q8, int M, int K);
105void gemv_q6_k_q8_k_avx(float *y, const void *W, const void *x_q8, int M, int K);
106void gemv_q6_k_q8_k_sse(float *y, const void *W, const void *x_q8, int M, int K);
107
108static int ck_q6k_q8k_force_ref(void)
109{
110 static int cached = -1;
111 if (cached < 0) {
112 const char *env = getenv("CK_DEBUG_Q6K_Q8K_REF");
113 cached = (env && env[0] && env[0] != '0') ? 1 : 0;
114 }
115 return cached;
116}
117
118/* ============================================================================
119 * Reference Implementation
120 * ============================================================================ */
121
122/**
123 * @brief Scalar dot product for Q6_K x Q8_K
124 *
125 * Q6_K layout: 256 weights per block
126 * - ql[0..127]: low 4 bits for all 256 weights (packed 2 per byte)
127 * - qh[0..63]: high 2 bits for all 256 weights (packed 4 per byte)
128 * - scales[0..15]: int8 scale for each 16-weight sub-block
129 * - d: FP16 super-block scale
130 *
131 * The dequantization formula for each weight is:
132 * weight = d * scale[sub] * (q6_value - 32)
133 * where q6_value is the 6-bit unsigned value (0..63).
134 */
135static float dot_q6_k_q8_k_ref(const block_q6_K *w,
136 const block_q8_K *x,
137 int K)
138{
139 const int nb = K / QK_K;
140 float sumf = 0.0f;
141 float sums[8] = {0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f, 0.0f};
142
143 for (int i = 0; i < nb; ++i) {
144 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
145
146 const uint8_t *ql = w[i].ql;
147 const uint8_t *qh = w[i].qh;
148 const int8_t *sc = w[i].scales;
149 const int8_t *q8 = x[i].qs;
150 int32_t aux32[8] = {0, 0, 0, 0, 0, 0, 0, 0};
151
152 /* Process 256 weights in 2 iterations of 128 */
153 for (int n = 0; n < QK_K; n += 128) {
154 /* Each iteration processes 128 weights:
155 * - ql[0..63] contains low 4 bits
156 * - qh[0..31] contains high 2 bits
157 * - Interleaved pattern: weights 0-31, 32-63, 64-95, 96-127
158 */
159 for (int l = 0; l < 32; ++l) {
160 /* Sub-block index: each scale covers 16 weights */
161 const int is = l / 16;
162
163 /* Extract 6-bit values from packed format */
164 /* q1: weights l+0 (low nibble of ql[l], bits 0-1 of qh[l]) */
165 const int8_t q1 = (int8_t)((ql[l + 0] & 0xF) | (((qh[l] >> 0) & 3) << 4)) - 32;
166 /* q2: weights l+32 (low nibble of ql[l+32], bits 2-3 of qh[l]) */
167 const int8_t q2 = (int8_t)((ql[l + 32] & 0xF) | (((qh[l] >> 2) & 3) << 4)) - 32;
168 /* q3: weights l+64 (high nibble of ql[l], bits 4-5 of qh[l]) */
169 const int8_t q3 = (int8_t)((ql[l + 0] >> 4) | (((qh[l] >> 4) & 3) << 4)) - 32;
170 /* q4: weights l+96 (high nibble of ql[l+32], bits 6-7 of qh[l]) */
171 const int8_t q4 = (int8_t)((ql[l + 32] >> 4) | (((qh[l] >> 6) & 3) << 4)) - 32;
172
173 aux32[l & 7] += (int)sc[is + 0] * (int)q1 * (int)q8[l + 0];
174 aux32[l & 7] += (int)sc[is + 2] * (int)q2 * (int)q8[l + 32];
175 aux32[l & 7] += (int)sc[is + 4] * (int)q3 * (int)q8[l + 64];
176 aux32[l & 7] += (int)sc[is + 6] * (int)q4 * (int)q8[l + 96];
177 }
178 q8 += 128;
179 ql += 64;
180 qh += 32;
181 sc += 8;
182 }
183
184 for (int l = 0; l < 8; ++l) {
185 sums[l] += d * (float)aux32[l];
186 }
187 }
188
189 for (int l = 0; l < 8; ++l) {
190 sumf += sums[l];
191 }
192 return sumf;
193}
194
195#if defined(__ARM_NEON) || defined(__aarch64__)
196static float dot_q6_k_q8_k_neon(const block_q6_K *w,
197 const block_q8_K *x,
198 int K)
199{
200 const int nb = K / QK_K;
201 float sumf = 0.0f;
202
203 for (int i = 0; i < nb; ++i) {
204 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
205
206 const uint8_t *ql = w[i].ql;
207 const uint8_t *qh = w[i].qh;
208 const int8_t *sc = w[i].scales;
209 const int8_t *q8 = x[i].qs;
210
211 int8_t wvals[QK_K];
212 int8_t svals[QK_K];
213
214 for (int n = 0; n < QK_K; n += 128) {
215 for (int l = 0; l < 32; ++l) {
216 const int is = l / 16;
217
218 const int8_t q1 = (int8_t)((ql[l + 0] & 0xF) | (((qh[l] >> 0) & 3) << 4)) - 32;
219 const int8_t q2 = (int8_t)((ql[l + 32] & 0xF) | (((qh[l] >> 2) & 3) << 4)) - 32;
220 const int8_t q3 = (int8_t)((ql[l + 0] >> 4) | (((qh[l] >> 4) & 3) << 4)) - 32;
221 const int8_t q4 = (int8_t)((ql[l + 32] >> 4) | (((qh[l] >> 6) & 3) << 4)) - 32;
222
223 const int base = n;
224 wvals[base + l + 0] = q1;
225 wvals[base + l + 32] = q2;
226 wvals[base + l + 64] = q3;
227 wvals[base + l + 96] = q4;
228
229 svals[base + l + 0] = sc[is + 0];
230 svals[base + l + 32] = sc[is + 2];
231 svals[base + l + 64] = sc[is + 4];
232 svals[base + l + 96] = sc[is + 6];
233 }
234
235 ql += 64;
236 qh += 32;
237 sc += 8;
238 }
239
240 int32x4_t acc = vdupq_n_s32(0);
241 for (int j = 0; j < QK_K; j += 16) {
242 const int8x16_t wv = vld1q_s8(&wvals[j]);
243 const int8x16_t sv = vld1q_s8(&svals[j]);
244 const int8x16_t xv = vld1q_s8(&q8[j]);
245
246 const int16x8_t ws0 = vmull_s8(vget_low_s8(wv), vget_low_s8(sv));
247 const int16x8_t ws1 = vmull_s8(vget_high_s8(wv), vget_high_s8(sv));
248 const int16x8_t x0 = vmovl_s8(vget_low_s8(xv));
249 const int16x8_t x1 = vmovl_s8(vget_high_s8(xv));
250
251 int32x4_t p0 = vmull_s16(vget_low_s16(ws0), vget_low_s16(x0));
252 p0 = vmlal_s16(p0, vget_high_s16(ws0), vget_high_s16(x0));
253
254 int32x4_t p1 = vmull_s16(vget_low_s16(ws1), vget_low_s16(x1));
255 p1 = vmlal_s16(p1, vget_high_s16(ws1), vget_high_s16(x1));
256
257 acc = vaddq_s32(acc, p0);
258 acc = vaddq_s32(acc, p1);
259 }
260
261 int32_t lanes[4];
262 vst1q_s32(lanes, acc);
263 sumf += d * (float)(lanes[0] + lanes[1] + lanes[2] + lanes[3]);
264 }
265
266 return sumf;
267}
268
269static void gemv_q6_k_q8_k_neon(float *y,
270 const void *W,
271 const void *x_q8,
272 int M, int K)
273{
274 const block_q6_K *blocks = (const block_q6_K *)W;
275 const block_q8_K *x = (const block_q8_K *)x_q8;
276 const int blocks_per_row = K / QK_K;
277
278 for (int row = 0; row < M; ++row) {
279 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
280 y[row] = dot_q6_k_q8_k_neon(w_row, x, K);
281 }
282}
283#endif
284
285void gemv_q6_k_q8_k_ref(float *y,
286 const void *W,
287 const void *x_q8,
288 int M, int K)
289{
290 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
291 return;
292 }
293
294 const block_q6_K *blocks = (const block_q6_K *)W;
295 const block_q8_K *x = (const block_q8_K *)x_q8;
296 const int blocks_per_row = K / QK_K;
297
298 for (int row = 0; row < M; ++row) {
299 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
300 y[row] = dot_q6_k_q8_k_ref(w_row, x, K);
301 }
302}
303
304/* ============================================================================
305 * SSE4.1 Implementation (for Ivy Bridge and older AVX-without-AVX2 CPUs)
306 *
307 * Uses 128-bit SSE operations with maddubs for integer multiply-add.
308 * Handles the -32 offset using bsums from Q8_K.
309 * ============================================================================ */
310
311#if defined(__SSSE3__)
312
313/* Scale shuffle indices for Q6_K - maps scale index to 16-byte shuffle pattern */
314static const int8_t q6k_scale_shuffle[8][16] = {
315 { 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1 }, /* is=0: scales[0,1] */
316 { 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 }, /* is=1: scales[2,3] */
317 { 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5 }, /* is=2: scales[4,5] */
318 { 6, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7 }, /* is=3: scales[6,7] */
319 { 8, 8, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 9 }, /* is=4: scales[8,9] */
320 {10,10,10,10,10,10,10,10,11,11,11,11,11,11,11,11 }, /* is=5: scales[10,11] */
321 {12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,13 }, /* is=6: scales[12,13] */
322 {14,14,14,14,14,14,14,14,15,15,15,15,15,15,15,15 }, /* is=7: scales[14,15] */
323};
324
325static float dot_q6_k_q8_k_sse(const block_q6_K *w,
326 const block_q8_K *x,
327 int K)
328{
329 const int nb = K / QK_K;
330 const __m128i m3 = _mm_set1_epi8(3);
331 const __m128i m15 = _mm_set1_epi8(15);
332
333 __m128 acc = _mm_setzero_ps();
334
335 for (int i = 0; i < nb; ++i) {
336 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
337
338 const uint8_t *ql = w[i].ql;
339 const uint8_t *qh = w[i].qh;
340 const int8_t *q8 = x[i].qs;
341
342 /* Load scales and precompute the -32 offset contribution using bsums */
343 const __m128i scales = _mm_loadu_si128((const __m128i *)w[i].scales);
344 const __m128i q8sums_0 = _mm_loadu_si128((const __m128i *)x[i].bsums);
345 const __m128i q8sums_1 = _mm_loadu_si128((const __m128i *)x[i].bsums + 1);
346
347 /* Compute: sum(scale * bsum) * 32 for the -32 offset */
348 const __m128i scales_16_0 = _mm_cvtepi8_epi16(scales);
349 const __m128i scales_16_1 = _mm_cvtepi8_epi16(_mm_bsrli_si128(scales, 8));
350 const __m128i q8sclsub_0 = _mm_slli_epi32(_mm_madd_epi16(q8sums_0, scales_16_0), 5);
351 const __m128i q8sclsub_1 = _mm_slli_epi32(_mm_madd_epi16(q8sums_1, scales_16_1), 5);
352
353 __m128i sumi_0 = _mm_setzero_si128();
354 __m128i sumi_1 = _mm_setzero_si128();
355
356 int is = 0;
357
358 /* Process 256 weights in 2 iterations of 128 */
359 for (int j = 0; j < QK_K / 128; ++j) {
360 /* Load high bits */
361 const __m128i q4bitsH_0 = _mm_loadu_si128((const __m128i *)qh);
362 qh += 16;
363 const __m128i q4bitsH_1 = _mm_loadu_si128((const __m128i *)qh);
364 qh += 16;
365
366 /* Extract and shift high bits into position */
367 const __m128i q4h_0 = _mm_slli_epi16(_mm_and_si128(q4bitsH_0, m3), 4);
368 const __m128i q4h_1 = _mm_slli_epi16(_mm_and_si128(q4bitsH_1, m3), 4);
369 const __m128i q4h_2 = _mm_slli_epi16(_mm_and_si128(q4bitsH_0, _mm_set1_epi8(12)), 2);
370 const __m128i q4h_3 = _mm_slli_epi16(_mm_and_si128(q4bitsH_1, _mm_set1_epi8(12)), 2);
371 const __m128i q4h_4 = _mm_and_si128(q4bitsH_0, _mm_set1_epi8(48));
372 const __m128i q4h_5 = _mm_and_si128(q4bitsH_1, _mm_set1_epi8(48));
373 const __m128i q4h_6 = _mm_srli_epi16(_mm_and_si128(q4bitsH_0, _mm_set1_epi8(-64)), 2);
374 const __m128i q4h_7 = _mm_srli_epi16(_mm_and_si128(q4bitsH_1, _mm_set1_epi8(-64)), 2);
375
376 /* Load low bits */
377 const __m128i q4bits1_0 = _mm_loadu_si128((const __m128i *)ql);
378 ql += 16;
379 const __m128i q4bits1_1 = _mm_loadu_si128((const __m128i *)ql);
380 ql += 16;
381 const __m128i q4bits2_0 = _mm_loadu_si128((const __m128i *)ql);
382 ql += 16;
383 const __m128i q4bits2_1 = _mm_loadu_si128((const __m128i *)ql);
384 ql += 16;
385
386 /* Combine low and high bits to get 6-bit values (unsigned 0..63) */
387 const __m128i q4_0 = _mm_or_si128(_mm_and_si128(q4bits1_0, m15), q4h_0);
388 const __m128i q4_1 = _mm_or_si128(_mm_and_si128(q4bits1_1, m15), q4h_1);
389 const __m128i q4_2 = _mm_or_si128(_mm_and_si128(q4bits2_0, m15), q4h_2);
390 const __m128i q4_3 = _mm_or_si128(_mm_and_si128(q4bits2_1, m15), q4h_3);
391 const __m128i q4_4 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits1_0, 4), m15), q4h_4);
392 const __m128i q4_5 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits1_1, 4), m15), q4h_5);
393 const __m128i q4_6 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits2_0, 4), m15), q4h_6);
394 const __m128i q4_7 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits2_1, 4), m15), q4h_7);
395
396 /* Load Q8_K values */
397 const __m128i q8_0 = _mm_loadu_si128((const __m128i *)q8);
398 q8 += 16;
399 const __m128i q8_1 = _mm_loadu_si128((const __m128i *)q8);
400 q8 += 16;
401 const __m128i q8_2 = _mm_loadu_si128((const __m128i *)q8);
402 q8 += 16;
403 const __m128i q8_3 = _mm_loadu_si128((const __m128i *)q8);
404 q8 += 16;
405 const __m128i q8_4 = _mm_loadu_si128((const __m128i *)q8);
406 q8 += 16;
407 const __m128i q8_5 = _mm_loadu_si128((const __m128i *)q8);
408 q8 += 16;
409 const __m128i q8_6 = _mm_loadu_si128((const __m128i *)q8);
410 q8 += 16;
411 const __m128i q8_7 = _mm_loadu_si128((const __m128i *)q8);
412 q8 += 16;
413
414 /* Multiply: maddubs treats first arg as unsigned, second as signed */
415 __m128i p16_0 = _mm_maddubs_epi16(q4_0, q8_0);
416 __m128i p16_1 = _mm_maddubs_epi16(q4_1, q8_1);
417 __m128i p16_2 = _mm_maddubs_epi16(q4_2, q8_2);
418 __m128i p16_3 = _mm_maddubs_epi16(q4_3, q8_3);
419 __m128i p16_4 = _mm_maddubs_epi16(q4_4, q8_4);
420 __m128i p16_5 = _mm_maddubs_epi16(q4_5, q8_5);
421 __m128i p16_6 = _mm_maddubs_epi16(q4_6, q8_6);
422 __m128i p16_7 = _mm_maddubs_epi16(q4_7, q8_7);
423
424 /* Get scales for this iteration */
425 const __m128i scale_0 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 0]));
426 const __m128i scale_1 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 1]));
427 const __m128i scale_2 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 2]));
428 const __m128i scale_3 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 3]));
429 is += 4;
430
431 /* Scale the products and widen to 32-bit */
432 p16_0 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_0), p16_0);
433 p16_1 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_0, 8)), p16_1);
434 p16_2 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_1), p16_2);
435 p16_3 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_1, 8)), p16_3);
436 p16_4 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_2), p16_4);
437 p16_5 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_2, 8)), p16_5);
438 p16_6 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_3), p16_6);
439 p16_7 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_3, 8)), p16_7);
440
441 /* Accumulate */
442 sumi_0 = _mm_add_epi32(sumi_0, _mm_add_epi32(p16_0, p16_2));
443 sumi_1 = _mm_add_epi32(sumi_1, _mm_add_epi32(p16_1, p16_3));
444 sumi_0 = _mm_add_epi32(sumi_0, _mm_add_epi32(p16_4, p16_6));
445 sumi_1 = _mm_add_epi32(sumi_1, _mm_add_epi32(p16_5, p16_7));
446 }
447
448 /* Subtract the -32 offset contribution */
449 sumi_0 = _mm_sub_epi32(sumi_0, q8sclsub_0);
450 sumi_1 = _mm_sub_epi32(sumi_1, q8sclsub_1);
451
452 /* Combine and convert to float */
453 __m128i sumi = _mm_add_epi32(sumi_0, sumi_1);
454 __m128 sumf_vec = _mm_mul_ps(_mm_set1_ps(d), _mm_cvtepi32_ps(sumi));
455
456 /* Horizontal sum */
457 sumf_vec = _mm_hadd_ps(sumf_vec, sumf_vec);
458 sumf_vec = _mm_hadd_ps(sumf_vec, sumf_vec);
459 acc = _mm_add_ss(acc, sumf_vec);
460 }
461
462 return _mm_cvtss_f32(acc);
463}
464
465void gemv_q6_k_q8_k_sse(float *y,
466 const void *W,
467 const void *x_q8,
468 int M, int K)
469{
470 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
471 return;
472 }
473
474 const block_q6_K *blocks = (const block_q6_K *)W;
475 const block_q8_K *x = (const block_q8_K *)x_q8;
476 const int blocks_per_row = K / QK_K;
477
478 for (int row = 0; row < M; ++row) {
479 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
480 y[row] = dot_q6_k_q8_k_sse(w_row, x, K);
481 }
482}
483#endif /* __SSSE3__ */
484
485/* ============================================================================
486 * AVX Implementation (for Sandy/Ivy Bridge - AVX without AVX2)
487 *
488 * Same as SSE but with prefetching for next block.
489 * Uses 128-bit integer ops (AVX doesn't add 256-bit int ops).
490 * ============================================================================ */
491
492#if defined(__AVX__) && !defined(__AVX2__)
493
494static float dot_q6_k_q8_k_avx(const block_q6_K *w,
495 const block_q8_K *x,
496 int K)
497{
498 const int nb = K / QK_K;
499 const __m128i m3 = _mm_set1_epi8(3);
500 const __m128i m15 = _mm_set1_epi8(15);
501
502 __m128 acc = _mm_setzero_ps();
503
504 for (int i = 0; i < nb; ++i) {
505 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
506
507 const uint8_t *ql = w[i].ql;
508 const uint8_t *qh = w[i].qh;
509 const int8_t *q8 = x[i].qs;
510
511 /* Prefetch next block */
512 if (i + 1 < nb) {
513 _mm_prefetch((const char *)&w[i + 1], _MM_HINT_T0);
514 _mm_prefetch((const char *)&x[i + 1], _MM_HINT_T0);
515 }
516
517 /* Load scales and precompute the -32 offset contribution using bsums */
518 const __m128i scales = _mm_loadu_si128((const __m128i *)w[i].scales);
519 const __m128i q8sums_0 = _mm_loadu_si128((const __m128i *)x[i].bsums);
520 const __m128i q8sums_1 = _mm_loadu_si128((const __m128i *)x[i].bsums + 1);
521
522 /* Compute: sum(scale * bsum) * 32 for the -32 offset */
523 const __m128i scales_16_0 = _mm_cvtepi8_epi16(scales);
524 const __m128i scales_16_1 = _mm_cvtepi8_epi16(_mm_bsrli_si128(scales, 8));
525 const __m128i q8sclsub_0 = _mm_slli_epi32(_mm_madd_epi16(q8sums_0, scales_16_0), 5);
526 const __m128i q8sclsub_1 = _mm_slli_epi32(_mm_madd_epi16(q8sums_1, scales_16_1), 5);
527
528 __m128i sumi_0 = _mm_setzero_si128();
529 __m128i sumi_1 = _mm_setzero_si128();
530
531 int is = 0;
532
533 /* Process 256 weights in 2 iterations of 128 */
534 for (int j = 0; j < QK_K / 128; ++j) {
535 /* Load high bits */
536 const __m128i q4bitsH_0 = _mm_loadu_si128((const __m128i *)qh);
537 qh += 16;
538 const __m128i q4bitsH_1 = _mm_loadu_si128((const __m128i *)qh);
539 qh += 16;
540
541 /* Extract and shift high bits into position */
542 const __m128i q4h_0 = _mm_slli_epi16(_mm_and_si128(q4bitsH_0, m3), 4);
543 const __m128i q4h_1 = _mm_slli_epi16(_mm_and_si128(q4bitsH_1, m3), 4);
544 const __m128i q4h_2 = _mm_slli_epi16(_mm_and_si128(q4bitsH_0, _mm_set1_epi8(12)), 2);
545 const __m128i q4h_3 = _mm_slli_epi16(_mm_and_si128(q4bitsH_1, _mm_set1_epi8(12)), 2);
546 const __m128i q4h_4 = _mm_and_si128(q4bitsH_0, _mm_set1_epi8(48));
547 const __m128i q4h_5 = _mm_and_si128(q4bitsH_1, _mm_set1_epi8(48));
548 const __m128i q4h_6 = _mm_srli_epi16(_mm_and_si128(q4bitsH_0, _mm_set1_epi8(-64)), 2);
549 const __m128i q4h_7 = _mm_srli_epi16(_mm_and_si128(q4bitsH_1, _mm_set1_epi8(-64)), 2);
550
551 /* Load low bits */
552 const __m128i q4bits1_0 = _mm_loadu_si128((const __m128i *)ql);
553 ql += 16;
554 const __m128i q4bits1_1 = _mm_loadu_si128((const __m128i *)ql);
555 ql += 16;
556 const __m128i q4bits2_0 = _mm_loadu_si128((const __m128i *)ql);
557 ql += 16;
558 const __m128i q4bits2_1 = _mm_loadu_si128((const __m128i *)ql);
559 ql += 16;
560
561 /* Combine low and high bits to get 6-bit values (unsigned 0..63) */
562 const __m128i q4_0 = _mm_or_si128(_mm_and_si128(q4bits1_0, m15), q4h_0);
563 const __m128i q4_1 = _mm_or_si128(_mm_and_si128(q4bits1_1, m15), q4h_1);
564 const __m128i q4_2 = _mm_or_si128(_mm_and_si128(q4bits2_0, m15), q4h_2);
565 const __m128i q4_3 = _mm_or_si128(_mm_and_si128(q4bits2_1, m15), q4h_3);
566 const __m128i q4_4 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits1_0, 4), m15), q4h_4);
567 const __m128i q4_5 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits1_1, 4), m15), q4h_5);
568 const __m128i q4_6 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits2_0, 4), m15), q4h_6);
569 const __m128i q4_7 = _mm_or_si128(_mm_and_si128(_mm_srli_epi16(q4bits2_1, 4), m15), q4h_7);
570
571 /* Load Q8_K values */
572 const __m128i q8_0 = _mm_loadu_si128((const __m128i *)q8);
573 q8 += 16;
574 const __m128i q8_1 = _mm_loadu_si128((const __m128i *)q8);
575 q8 += 16;
576 const __m128i q8_2 = _mm_loadu_si128((const __m128i *)q8);
577 q8 += 16;
578 const __m128i q8_3 = _mm_loadu_si128((const __m128i *)q8);
579 q8 += 16;
580 const __m128i q8_4 = _mm_loadu_si128((const __m128i *)q8);
581 q8 += 16;
582 const __m128i q8_5 = _mm_loadu_si128((const __m128i *)q8);
583 q8 += 16;
584 const __m128i q8_6 = _mm_loadu_si128((const __m128i *)q8);
585 q8 += 16;
586 const __m128i q8_7 = _mm_loadu_si128((const __m128i *)q8);
587 q8 += 16;
588
589 /* Multiply: maddubs treats first arg as unsigned, second as signed */
590 __m128i p16_0 = _mm_maddubs_epi16(q4_0, q8_0);
591 __m128i p16_1 = _mm_maddubs_epi16(q4_1, q8_1);
592 __m128i p16_2 = _mm_maddubs_epi16(q4_2, q8_2);
593 __m128i p16_3 = _mm_maddubs_epi16(q4_3, q8_3);
594 __m128i p16_4 = _mm_maddubs_epi16(q4_4, q8_4);
595 __m128i p16_5 = _mm_maddubs_epi16(q4_5, q8_5);
596 __m128i p16_6 = _mm_maddubs_epi16(q4_6, q8_6);
597 __m128i p16_7 = _mm_maddubs_epi16(q4_7, q8_7);
598
599 /* Get scales for this iteration */
600 const __m128i scale_0 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 0]));
601 const __m128i scale_1 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 1]));
602 const __m128i scale_2 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 2]));
603 const __m128i scale_3 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)q6k_scale_shuffle[is + 3]));
604 is += 4;
605
606 /* Scale the products and widen to 32-bit */
607 p16_0 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_0), p16_0);
608 p16_1 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_0, 8)), p16_1);
609 p16_2 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_1), p16_2);
610 p16_3 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_1, 8)), p16_3);
611 p16_4 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_2), p16_4);
612 p16_5 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_2, 8)), p16_5);
613 p16_6 = _mm_madd_epi16(_mm_cvtepi8_epi16(scale_3), p16_6);
614 p16_7 = _mm_madd_epi16(_mm_cvtepi8_epi16(_mm_bsrli_si128(scale_3, 8)), p16_7);
615
616 /* Accumulate */
617 sumi_0 = _mm_add_epi32(sumi_0, _mm_add_epi32(p16_0, p16_2));
618 sumi_1 = _mm_add_epi32(sumi_1, _mm_add_epi32(p16_1, p16_3));
619 sumi_0 = _mm_add_epi32(sumi_0, _mm_add_epi32(p16_4, p16_6));
620 sumi_1 = _mm_add_epi32(sumi_1, _mm_add_epi32(p16_5, p16_7));
621 }
622
623 /* Subtract the -32 offset contribution */
624 sumi_0 = _mm_sub_epi32(sumi_0, q8sclsub_0);
625 sumi_1 = _mm_sub_epi32(sumi_1, q8sclsub_1);
626
627 /* Combine and convert to float */
628 __m128i sumi = _mm_add_epi32(sumi_0, sumi_1);
629 __m128 sumf_vec = _mm_mul_ps(_mm_set1_ps(d), _mm_cvtepi32_ps(sumi));
630
631 /* Horizontal sum */
632 sumf_vec = _mm_hadd_ps(sumf_vec, sumf_vec);
633 sumf_vec = _mm_hadd_ps(sumf_vec, sumf_vec);
634 acc = _mm_add_ss(acc, sumf_vec);
635 }
636
637 return _mm_cvtss_f32(acc);
638}
639
640void gemv_q6_k_q8_k_avx(float *y,
641 const void *W,
642 const void *x_q8,
643 int M, int K)
644{
645 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
646 return;
647 }
648
649 const block_q6_K *blocks = (const block_q6_K *)W;
650 const block_q8_K *x = (const block_q8_K *)x_q8;
651 const int blocks_per_row = K / QK_K;
652
653 for (int row = 0; row < M; ++row) {
654 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
655 y[row] = dot_q6_k_q8_k_avx(w_row, x, K);
656 }
657}
658
659#endif /* __AVX__ && !__AVX2__ */
660
661/* ============================================================================
662 * AVX2 Implementation (for modern CPUs with AVX2)
663 * ============================================================================ */
664
665#if defined(__AVX2__)
666
667/* Match ggml's x86 hsum_float_8 reduction tree exactly. Replacing the final
668 * two additions with _mm_hadd_ps changes one FP32 ULP for realistic MLP-down
669 * rows and can cross the absolute parity gate on AVX2/AVX-512 runners. */
670static inline float ck_q6k_hsum_float_8(const __m256 x)
671{
672 __m128 res = _mm256_extractf128_ps(x, 1);
673 res = _mm_add_ps(res, _mm256_castps256_ps128(x));
674 res = _mm_add_ps(res, _mm_movehl_ps(res, res));
675 res = _mm_add_ss(res, _mm_movehdup_ps(res));
676 return _mm_cvtss_f32(res);
677}
678
679/* Q6_K optimization note:
680 * ----------------------
681 * llama.cpp's newer Q6_K AVX2 work inspired an experiment that separates the
682 * zero-point correction from the hot unpack/dot loop by using Q8_K block sums
683 * (bsums). That can reduce repeated subtract-32 work and improve some CPU
684 * layouts, but it also changes reduction order and instruction balance.
685 *
686 * Current llama.cpp x86 kernels use the AVX2 reduction tree even in AVX-512
687 * builds. CK uses the same tree for its AVX2 and non-VBMI AVX-512 paths. Keep
688 * the final horizontal reduction in sync with ggml: an equivalent _mm_hadd_ps
689 * tree differs by one FP32 ULP on realistic MLP-down rows.
690 */
691
692static inline __m128i get_scale_shuffle_avx2(int i) {
693 static const uint8_t patterns[8][16] = {
694 { 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1 },
695 { 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 },
696 { 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5 },
697 { 6, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7 },
698 { 8, 8, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 9 },
699 {10,10,10,10,10,10,10,10,11,11,11,11,11,11,11,11 },
700 {12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,13 },
701 {14,14,14,14,14,14,14,14,15,15,15,15,15,15,15,15 },
702 };
703 return _mm_loadu_si128((const __m128i *)patterns[i]);
704}
705
706static float dot_q6_k_q8_k_avx2(const block_q6_K *w,
707 const block_q8_K *x,
708 int K)
709{
710 const int nb = K / QK_K;
711 const __m256i m3 = _mm256_set1_epi8(3);
712 const __m256i m15 = _mm256_set1_epi8(15);
713
714 __m256 acc = _mm256_setzero_ps();
715
716 for (int i = 0; i < nb; ++i) {
717 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
718
719 const uint8_t *q4 = w[i].ql;
720 const uint8_t *qh = w[i].qh;
721 const int8_t *q8 = x[i].qs;
722
723 const __m256i q8sums = _mm256_loadu_si256((const __m256i *)x[i].bsums);
724 const __m128i scales = _mm_loadu_si128((const __m128i *)w[i].scales);
725 const __m256i scales_16 = _mm256_cvtepi8_epi16(scales);
726 const __m256i q8sclsub = _mm256_slli_epi32(
727 _mm256_madd_epi16(q8sums, scales_16), 5);
728
729 __m256i sumi = _mm256_setzero_si256();
730 int is = 0;
731
732 for (int j = 0; j < QK_K / 128; ++j) {
733 const __m256i q4bits1 = _mm256_loadu_si256((const __m256i *)q4);
734 q4 += 32;
735 const __m256i q4bits2 = _mm256_loadu_si256((const __m256i *)q4);
736 q4 += 32;
737 const __m256i q4bitsH = _mm256_loadu_si256((const __m256i *)qh);
738 qh += 32;
739
740 const __m256i q4h_0 = _mm256_slli_epi16(_mm256_and_si256(q4bitsH, m3), 4);
741 const __m256i q4h_1 = _mm256_slli_epi16(
742 _mm256_and_si256(q4bitsH, _mm256_set1_epi8(12)), 2);
743 const __m256i q4h_2 = _mm256_and_si256(q4bitsH, _mm256_set1_epi8(48));
744 const __m256i q4h_3 = _mm256_srli_epi16(
745 _mm256_and_si256(q4bitsH, _mm256_set1_epi8(-64)), 2);
746
747 const __m256i q4_0 = _mm256_or_si256(_mm256_and_si256(q4bits1, m15), q4h_0);
748 const __m256i q4_1 = _mm256_or_si256(_mm256_and_si256(q4bits2, m15), q4h_1);
749 const __m256i q4_2 = _mm256_or_si256(
750 _mm256_and_si256(_mm256_srli_epi16(q4bits1, 4), m15), q4h_2);
751 const __m256i q4_3 = _mm256_or_si256(
752 _mm256_and_si256(_mm256_srli_epi16(q4bits2, 4), m15), q4h_3);
753
754 const __m256i q8_0 = _mm256_loadu_si256((const __m256i *)q8);
755 q8 += 32;
756 const __m256i q8_1 = _mm256_loadu_si256((const __m256i *)q8);
757 q8 += 32;
758 const __m256i q8_2 = _mm256_loadu_si256((const __m256i *)q8);
759 q8 += 32;
760 const __m256i q8_3 = _mm256_loadu_si256((const __m256i *)q8);
761 q8 += 32;
762
763 __m256i p16_0 = _mm256_maddubs_epi16(q4_0, q8_0);
764 __m256i p16_1 = _mm256_maddubs_epi16(q4_1, q8_1);
765 __m256i p16_2 = _mm256_maddubs_epi16(q4_2, q8_2);
766 __m256i p16_3 = _mm256_maddubs_epi16(q4_3, q8_3);
767
768 const __m128i scale_0 = _mm_shuffle_epi8(
769 scales, get_scale_shuffle_avx2(is + 0));
770 const __m128i scale_1 = _mm_shuffle_epi8(
771 scales, get_scale_shuffle_avx2(is + 1));
772 const __m128i scale_2 = _mm_shuffle_epi8(
773 scales, get_scale_shuffle_avx2(is + 2));
774 const __m128i scale_3 = _mm_shuffle_epi8(
775 scales, get_scale_shuffle_avx2(is + 3));
776 is += 4;
777
778 p16_0 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_0), p16_0);
779 p16_1 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_1), p16_1);
780 p16_2 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_2), p16_2);
781 p16_3 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_3), p16_3);
782
783 sumi = _mm256_add_epi32(sumi, _mm256_add_epi32(p16_0, p16_1));
784 sumi = _mm256_add_epi32(sumi, _mm256_add_epi32(p16_2, p16_3));
785 }
786
787 sumi = _mm256_sub_epi32(sumi, q8sclsub);
788 acc = _mm256_fmadd_ps(_mm256_broadcast_ss(&d), _mm256_cvtepi32_ps(sumi), acc);
789 }
790
791 return ck_q6k_hsum_float_8(acc);
792}
793
794/*
795 * Four-token Q6_K x Q8_K microkernel.
796 *
797 * Each output keeps the same independent integer and FP32 reduction tree as
798 * dot_q6_k_q8_k_avx2(). The only changed ordering is between independent
799 * outputs: Q6 low/high-bit unpack and scale shuffles are performed once, then
800 * reused by as many as four Q8_K activation rows. This makes the provider
801 * eligible for exact llama.cpp parity while removing repeated Q6 decode work.
802 */
803static void dot_q6_k_q8_k_avx2_m4(const block_q6_K *w,
804 const block_q8_K *x,
805 int x_row_blocks,
806 int rows,
807 int K,
808 float out[4])
809{
810 const int nb = K / QK_K;
811 const __m256i m3 = _mm256_set1_epi8(3);
812 const __m256i m15 = _mm256_set1_epi8(15);
813 __m256 acc[4] = {
814 _mm256_setzero_ps(), _mm256_setzero_ps(),
815 _mm256_setzero_ps(), _mm256_setzero_ps()
816 };
817
818 for (int i = 0; i < nb; ++i) {
819 const uint8_t *q4 = w[i].ql;
820 const uint8_t *qh = w[i].qh;
821 const __m128i scales = _mm_loadu_si128((const __m128i *)w[i].scales);
822 const __m256i scales_16 = _mm256_cvtepi8_epi16(scales);
823 const float wd = GGML_FP16_TO_FP32(w[i].d);
824 __m256i sumi[4] = {
825 _mm256_setzero_si256(), _mm256_setzero_si256(),
826 _mm256_setzero_si256(), _mm256_setzero_si256()
827 };
828 int is = 0;
829
830 for (int j = 0; j < QK_K / 128; ++j) {
831 const __m256i q4bits1 = _mm256_loadu_si256((const __m256i *)q4);
832 q4 += 32;
833 const __m256i q4bits2 = _mm256_loadu_si256((const __m256i *)q4);
834 q4 += 32;
835 const __m256i q4bitsH = _mm256_loadu_si256((const __m256i *)qh);
836 qh += 32;
837
838 const __m256i q4h_0 = _mm256_slli_epi16(_mm256_and_si256(q4bitsH, m3), 4);
839 const __m256i q4h_1 = _mm256_slli_epi16(
840 _mm256_and_si256(q4bitsH, _mm256_set1_epi8(12)), 2);
841 const __m256i q4h_2 = _mm256_and_si256(q4bitsH, _mm256_set1_epi8(48));
842 const __m256i q4h_3 = _mm256_srli_epi16(
843 _mm256_and_si256(q4bitsH, _mm256_set1_epi8(-64)), 2);
844
845 const __m256i q6_0 = _mm256_or_si256(_mm256_and_si256(q4bits1, m15), q4h_0);
846 const __m256i q6_1 = _mm256_or_si256(_mm256_and_si256(q4bits2, m15), q4h_1);
847 const __m256i q6_2 = _mm256_or_si256(
848 _mm256_and_si256(_mm256_srli_epi16(q4bits1, 4), m15), q4h_2);
849 const __m256i q6_3 = _mm256_or_si256(
850 _mm256_and_si256(_mm256_srli_epi16(q4bits2, 4), m15), q4h_3);
851
852 const __m128i scale_0 = _mm_shuffle_epi8(
853 scales, get_scale_shuffle_avx2(is + 0));
854 const __m128i scale_1 = _mm_shuffle_epi8(
855 scales, get_scale_shuffle_avx2(is + 1));
856 const __m128i scale_2 = _mm_shuffle_epi8(
857 scales, get_scale_shuffle_avx2(is + 2));
858 const __m128i scale_3 = _mm_shuffle_epi8(
859 scales, get_scale_shuffle_avx2(is + 3));
860 is += 4;
861 const __m256i scale16_0 = _mm256_cvtepi8_epi16(scale_0);
862 const __m256i scale16_1 = _mm256_cvtepi8_epi16(scale_1);
863 const __m256i scale16_2 = _mm256_cvtepi8_epi16(scale_2);
864 const __m256i scale16_3 = _mm256_cvtepi8_epi16(scale_3);
865
866 for (int r = 0; r < rows; ++r) {
867 const block_q8_K *xb = x + (size_t)r * (size_t)x_row_blocks + i;
868 const int8_t *q8 = xb->qs + j * 128;
869 const __m256i q8_0 = _mm256_loadu_si256((const __m256i *)(q8 + 0));
870 const __m256i q8_1 = _mm256_loadu_si256((const __m256i *)(q8 + 32));
871 const __m256i q8_2 = _mm256_loadu_si256((const __m256i *)(q8 + 64));
872 const __m256i q8_3 = _mm256_loadu_si256((const __m256i *)(q8 + 96));
873
874 __m256i p16_0 = _mm256_maddubs_epi16(q6_0, q8_0);
875 __m256i p16_1 = _mm256_maddubs_epi16(q6_1, q8_1);
876 __m256i p16_2 = _mm256_maddubs_epi16(q6_2, q8_2);
877 __m256i p16_3 = _mm256_maddubs_epi16(q6_3, q8_3);
878 p16_0 = _mm256_madd_epi16(scale16_0, p16_0);
879 p16_1 = _mm256_madd_epi16(scale16_1, p16_1);
880 p16_2 = _mm256_madd_epi16(scale16_2, p16_2);
881 p16_3 = _mm256_madd_epi16(scale16_3, p16_3);
882
883 sumi[r] = _mm256_add_epi32(
884 sumi[r], _mm256_add_epi32(p16_0, p16_1));
885 sumi[r] = _mm256_add_epi32(
886 sumi[r], _mm256_add_epi32(p16_2, p16_3));
887 }
888 }
889
890 for (int r = 0; r < rows; ++r) {
891 const block_q8_K *xb = x + (size_t)r * (size_t)x_row_blocks + i;
892 const __m256i q8sums =
893 _mm256_loadu_si256((const __m256i *)xb->bsums);
894 const __m256i q8sclsub = _mm256_slli_epi32(
895 _mm256_madd_epi16(q8sums, scales_16), 5);
896 sumi[r] = _mm256_sub_epi32(sumi[r], q8sclsub);
897 const float d = wd * xb->d;
898 acc[r] = _mm256_fmadd_ps(
899 _mm256_broadcast_ss(&d), _mm256_cvtepi32_ps(sumi[r]), acc[r]);
900 }
901 }
902
903 for (int r = 0; r < rows; ++r) {
904 out[r] = ck_q6k_hsum_float_8(acc[r]);
905 }
906}
907
908static float dot_q6_k_prepared_q8_k_avx2(
909 const block_q6_K_prepared *w, const block_q8_K *x, int K)
910{
911 const int nb = K / QK_K;
912 __m256 acc = _mm256_setzero_ps();
913
914 for (int i = 0; i < nb; ++i) {
915 const __m128i scales =
916 _mm_loadu_si128((const __m128i *)(const void *)w[i].scales);
917 const __m256i scales_16 = _mm256_cvtepi8_epi16(scales);
918 const __m256i q8sums =
919 _mm256_loadu_si256((const __m256i *)(const void *)x[i].bsums);
920 const __m256i q8sclsub = _mm256_slli_epi32(
921 _mm256_madd_epi16(q8sums, scales_16), 5);
922 __m256i sumi = _mm256_setzero_si256();
923
924 for (int j = 0; j < 2; ++j) {
925 __m256i p16[4];
926 for (int lane = 0; lane < 4; ++lane) {
927 const int group = j * 4 + lane;
928 const __m256i q6 = _mm256_loadu_si256(
929 (const __m256i *)(const void *)(w[i].qs + group * 32));
930 const __m256i q8 = _mm256_loadu_si256(
931 (const __m256i *)(const void *)(x[i].qs + group * 32));
932 p16[lane] = _mm256_maddubs_epi16(q6, q8);
933 const __m128i scale = _mm_shuffle_epi8(
934 scales, get_scale_shuffle_avx2(group));
935 p16[lane] = _mm256_madd_epi16(
936 _mm256_cvtepi8_epi16(scale), p16[lane]);
937 }
938 sumi = _mm256_add_epi32(
939 sumi, _mm256_add_epi32(p16[0], p16[1]));
940 sumi = _mm256_add_epi32(
941 sumi, _mm256_add_epi32(p16[2], p16[3]));
942 }
943
944 sumi = _mm256_sub_epi32(sumi, q8sclsub);
945 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
946 acc = _mm256_fmadd_ps(
947 _mm256_broadcast_ss(&d), _mm256_cvtepi32_ps(sumi), acc);
948 }
949
950 return ck_q6k_hsum_float_8(acc);
951}
952
953#if defined(__AVX512F__) && defined(__AVX512BW__) && \
954 defined(__AVX512VNNI__)
955/*
956 * Prepared Q6 stores unsigned 6-bit values as bytes, which maps directly to
957 * VPDPBUSD. Fold each pair of 32-byte groups into the same eight integer lanes
958 * used by the certified AVX2 provider. The correction, FP32 block FMA, and
959 * horizontal reduction therefore retain the established arithmetic contract.
960 */
961static float dot_q6_k_prepared_q8_k_avx512_vnni(
962 const block_q6_K_prepared *w, const block_q8_K *x, int K)
963{
964 const int nb = K / QK_K;
965 const __m512i scale_indices = _mm512_setr_epi32(
966 0, 0, 0, 0, 1, 1, 1, 1,
967 2, 2, 2, 2, 3, 3, 3, 3);
968 __m256 acc = _mm256_setzero_ps();
969
970 for (int i = 0; i < nb; ++i) {
971 const __m128i scales =
972 _mm_loadu_si128((const __m128i *)(const void *)w[i].scales);
973 const __m512i scales_32 = _mm512_cvtepi8_epi32(scales);
974 const __m256i scales_16 = _mm256_cvtepi8_epi16(scales);
975 const __m256i q8sums =
976 _mm256_loadu_si256((const __m256i *)(const void *)x[i].bsums);
977 const __m256i q8sclsub = _mm256_slli_epi32(
978 _mm256_madd_epi16(q8sums, scales_16), 5);
979 __m256i sumi = _mm256_setzero_si256();
980
981 for (int chunk = 0; chunk < 4; ++chunk) {
982 const __m512i q6 = _mm512_loadu_si512(
983 (const void *)(w[i].qs + chunk * 64));
984 const __m512i q8 = _mm512_loadu_si512(
985 (const void *)(x[i].qs + chunk * 64));
986 const __m512i products = _mm512_dpbusd_epi32(
987 _mm512_setzero_si512(), q6, q8);
988 const __m512i indices = _mm512_add_epi32(
989 scale_indices, _mm512_set1_epi32(chunk * 4));
990 const __m512i scaled = _mm512_mullo_epi32(
991 products, _mm512_permutexvar_epi32(indices, scales_32));
992 const __m512i upper_half = _mm512_shuffle_i32x4(
993 scaled, scaled, _MM_SHUFFLE(3, 2, 3, 2));
994 const __m256i pair_sum = _mm256_add_epi32(
995 _mm512_castsi512_si256(scaled),
996 _mm512_castsi512_si256(upper_half));
997 sumi = _mm256_add_epi32(sumi, pair_sum);
998 }
999
1000 sumi = _mm256_sub_epi32(sumi, q8sclsub);
1001 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
1002 acc = _mm256_fmadd_ps(
1003 _mm256_broadcast_ss(&d), _mm256_cvtepi32_ps(sumi), acc);
1004 }
1005
1006 return ck_q6k_hsum_float_8(acc);
1007}
1008#endif
1009
1010void gemv_q6_k_q8_k_avx2(float *y,
1011 const void *W,
1012 const void *x_q8,
1013 int M, int K)
1014{
1015 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
1016 return;
1017 }
1018
1019 const block_q6_K *blocks = (const block_q6_K *)W;
1020 const block_q8_K *x = (const block_q8_K *)x_q8;
1021 const int blocks_per_row = K / QK_K;
1022
1023 for (int row = 0; row < M; ++row) {
1024 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
1025 y[row] = dot_q6_k_q8_k_avx2(w_row, x, K);
1026 }
1027}
1028#endif /* __AVX2__ */
1029
1030/* ============================================================================
1031 * AVX-512 Implementation
1032 *
1033 * Uses 512-bit ZMM registers to process 64 bytes at a time.
1034 * Processes entire 256-element Q6_K block in fewer iterations.
1035 * ============================================================================ */
1036
1037#if defined(__AVX512F__) && defined(__AVX512BW__) && defined(__AVX512VBMI__)
1038
1039/**
1040 * @brief AVX-512 dot product for Q6_K x Q8_K with VBMI
1041 *
1042 * Uses AVX-512 VBMI for efficient byte permutation.
1043 */
1044static float dot_q6_k_q8_k_avx512_vbmi(const block_q6_K *w,
1045 const block_q8_K *x,
1046 int K)
1047{
1048 const int nb = K / QK_K;
1049 const __m512i m4 = _mm512_set1_epi8(0xF);
1050 const __m512i m2 = _mm512_set1_epi8(3);
1051 const __m512i m32s = _mm512_set1_epi8(32);
1052
1053 __m512 acc = _mm512_setzero_ps();
1054
1055 for (int i = 0; i < nb; ++i) {
1056 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
1057
1058 const uint8_t *ql = w[i].ql;
1059 const uint8_t *qh = w[i].qh;
1060 const int8_t *q8 = x[i].qs;
1061 const int8_t *sc = w[i].scales;
1062
1063 __m512i sumi = _mm512_setzero_si512();
1064
1065 /* Process 256 weights in one iteration using AVX-512 */
1066 /* Load 64 bytes of low bits (covers 128 weights, need 2 loads for full block) */
1067 const __m512i q4bits1 = _mm512_loadu_si512((const __m512i *)ql); /* ql[0..63] */
1068 const __m512i q4bits2 = _mm512_loadu_si512((const __m512i *)(ql + 64)); /* ql[64..127] */
1069
1070 /* Load 64 bytes of high bits */
1071 const __m512i q4bitsH = _mm512_loadu_si512((const __m512i *)qh);
1072
1073 /* Extract high 2-bit contributions for each group of 32 weights */
1074 /* Group 0: bits 0-1 of qh -> weights 0-31 */
1075 const __m512i q4h_0 = _mm512_slli_epi16(_mm512_and_si512(q4bitsH, m2), 4);
1076 /* Group 1: bits 2-3 of qh -> weights 32-63 */
1077 const __m512i q4h_1 = _mm512_slli_epi16(_mm512_and_si512(_mm512_srli_epi16(q4bitsH, 2), m2), 4);
1078 /* Group 2: bits 4-5 of qh -> weights 64-95 */
1079 const __m512i q4h_2 = _mm512_slli_epi16(_mm512_and_si512(_mm512_srli_epi16(q4bitsH, 4), m2), 4);
1080 /* Group 3: bits 6-7 of qh -> weights 96-127 */
1081 const __m512i q4h_3 = _mm512_slli_epi16(_mm512_and_si512(_mm512_srli_epi16(q4bitsH, 6), m2), 4);
1082
1083 /* Combine low nibbles with high bits to get 6-bit values (0-63) */
1084 /* First 64 weights: low nibbles of ql[0..63] */
1085 const __m512i q6_0 = _mm512_or_si512(_mm512_and_si512(q4bits1, m4), q4h_0);
1086 const __m512i q6_1 = _mm512_or_si512(_mm512_and_si512(q4bits2, m4), q4h_1);
1087 /* Second 64 weights: high nibbles of ql[0..63] */
1088 const __m512i q6_2 = _mm512_or_si512(_mm512_and_si512(_mm512_srli_epi16(q4bits1, 4), m4), q4h_2);
1089 const __m512i q6_3 = _mm512_or_si512(_mm512_and_si512(_mm512_srli_epi16(q4bits2, 4), m4), q4h_3);
1090
1091 /* Load Q8_K values (256 int8 values = 4 x 64) */
1092 const __m512i q8_0 = _mm512_loadu_si512((const __m512i *)q8);
1093 const __m512i q8_1 = _mm512_loadu_si512((const __m512i *)(q8 + 64));
1094 const __m512i q8_2 = _mm512_loadu_si512((const __m512i *)(q8 + 128));
1095 const __m512i q8_3 = _mm512_loadu_si512((const __m512i *)(q8 + 192));
1096
1097 /* Compute 32 * q8 for the offset subtraction */
1098 __m512i q8s_0 = _mm512_maddubs_epi16(m32s, q8_0);
1099 __m512i q8s_1 = _mm512_maddubs_epi16(m32s, q8_1);
1100 __m512i q8s_2 = _mm512_maddubs_epi16(m32s, q8_2);
1101 __m512i q8s_3 = _mm512_maddubs_epi16(m32s, q8_3);
1102
1103 /* Multiply unsigned q6 * signed q8 */
1104 __m512i p16_0 = _mm512_maddubs_epi16(q6_0, q8_0);
1105 __m512i p16_1 = _mm512_maddubs_epi16(q6_1, q8_1);
1106 __m512i p16_2 = _mm512_maddubs_epi16(q6_2, q8_2);
1107 __m512i p16_3 = _mm512_maddubs_epi16(q6_3, q8_3);
1108
1109 /* Subtract offset: (q6 - 32) * q8 = q6*q8 - 32*q8 */
1110 p16_0 = _mm512_sub_epi16(p16_0, q8s_0);
1111 p16_1 = _mm512_sub_epi16(p16_1, q8s_1);
1112 p16_2 = _mm512_sub_epi16(p16_2, q8s_2);
1113 p16_3 = _mm512_sub_epi16(p16_3, q8s_3);
1114
1115 /* Load and broadcast scales using VBMI permute
1116 * Each scale applies to 16 weights, so we need to broadcast appropriately
1117 * scales[0..15] for the 16 sub-blocks */
1118 const __m128i scales_128 = _mm_loadu_si128((const __m128i *)sc);
1119
1120 /* Create scale broadcast patterns for 64 weights (4 scales per 64 weights) */
1121 /* Pattern: each scale repeated 16 times for 16 weights */
1122 const __m512i scale_idx_0 = _mm512_set_epi8(
1123 3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,
1124 2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,2,
1125 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,
1126 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0);
1127 const __m512i scale_idx_1 = _mm512_set_epi8(
1128 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
1129 6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
1130 5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
1131 4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4);
1132 const __m512i scale_idx_2 = _mm512_set_epi8(
1133 11,11,11,11,11,11,11,11,11,11,11,11,11,11,11,11,
1134 10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,10,
1135 9,9,9,9,9,9,9,9,9,9,9,9,9,9,9,9,
1136 8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8);
1137 const __m512i scale_idx_3 = _mm512_set_epi8(
1138 15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,15,
1139 14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,14,
1140 13,13,13,13,13,13,13,13,13,13,13,13,13,13,13,13,
1141 12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12);
1142
1143 /* Broadcast scales to 512-bit using VBMI permutexvar */
1144 const __m512i scales_512 = _mm512_broadcast_i32x4(scales_128);
1145 const __m512i sc_0 = _mm512_permutexvar_epi8(scale_idx_0, scales_512);
1146 const __m512i sc_1 = _mm512_permutexvar_epi8(scale_idx_1, scales_512);
1147 const __m512i sc_2 = _mm512_permutexvar_epi8(scale_idx_2, scales_512);
1148 const __m512i sc_3 = _mm512_permutexvar_epi8(scale_idx_3, scales_512);
1149
1150 /* Sign-extend scales to 16-bit and multiply with products */
1151 /* For efficiency, we process in two halves (low and high 256 bits) */
1152 __m512i p32_0 = _mm512_madd_epi16(_mm512_cvtepi8_epi16(_mm512_castsi512_si256(sc_0)), p16_0);
1153 __m512i p32_1 = _mm512_madd_epi16(_mm512_cvtepi8_epi16(_mm512_castsi512_si256(sc_1)), p16_1);
1154 __m512i p32_2 = _mm512_madd_epi16(_mm512_cvtepi8_epi16(_mm512_castsi512_si256(sc_2)), p16_2);
1155 __m512i p32_3 = _mm512_madd_epi16(_mm512_cvtepi8_epi16(_mm512_castsi512_si256(sc_3)), p16_3);
1156
1157 /* Accumulate */
1158 sumi = _mm512_add_epi32(sumi, p32_0);
1159 sumi = _mm512_add_epi32(sumi, p32_1);
1160 sumi = _mm512_add_epi32(sumi, p32_2);
1161 sumi = _mm512_add_epi32(sumi, p32_3);
1162
1163 /* Scale by d and accumulate */
1164 acc = _mm512_fmadd_ps(_mm512_set1_ps(d), _mm512_cvtepi32_ps(sumi), acc);
1165 }
1166
1167 return _mm512_reduce_add_ps(acc);
1168}
1169
1170void gemv_q6_k_q8_k_avx512_vbmi(float *y,
1171 const void *W,
1172 const void *x_q8,
1173 int M, int K)
1174{
1175 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
1176 return;
1177 }
1178
1179 const block_q6_K *blocks = (const block_q6_K *)W;
1180 const block_q8_K *x = (const block_q8_K *)x_q8;
1181 const int blocks_per_row = K / QK_K;
1182
1183 for (int row = 0; row < M; ++row) {
1184 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
1185 y[row] = dot_q6_k_q8_k_avx512_vbmi(w_row, x, K);
1186 }
1187}
1188
1189#endif /* __AVX512F__ && __AVX512BW__ && __AVX512VBMI__ */
1190
1191#if defined(__AVX512F__) && defined(__AVX512BW__)
1192
1193/**
1194 * @brief AVX-512 dot product for Q6_K x Q8_K
1195 *
1196 * Works on all AVX-512 CPUs (Skylake-X and newer).
1197 * Uses same algorithm as AVX2, but benefits from AVX-512's wider FMA
1198 * and efficient horizontal reduction.
1199 */
1200static float dot_q6_k_q8_k_avx512(const block_q6_K *w,
1201 const block_q8_K *x,
1202 int K)
1203{
1204 const int nb = K / QK_K;
1205 const __m256i m4 = _mm256_set1_epi8(0xF);
1206 const __m256i m2 = _mm256_set1_epi8(3);
1207 const __m256i m32s = _mm256_set1_epi8(32);
1208
1209 /* Use 256-bit float accumulator, same as AVX2 */
1210 __m256 acc = _mm256_setzero_ps();
1211
1212 for (int i = 0; i < nb; ++i) {
1213 const float d = GGML_FP16_TO_FP32(w[i].d) * x[i].d;
1214
1215 const uint8_t *q4 = w[i].ql;
1216 const uint8_t *qh = w[i].qh;
1217 const int8_t *q8 = x[i].qs;
1218
1219 const __m128i scales = _mm_loadu_si128((const __m128i *)w[i].scales);
1220
1221 /* Use 256-bit integer accumulator, same as AVX2 */
1222 __m256i sumi = _mm256_setzero_si256();
1223 int is = 0;
1224
1225 /* Process 256 weights in 2 iterations of 128 (same structure as AVX2) */
1226 for (int j = 0; j < QK_K / 128; ++j) {
1227 /* Get scale shuffle patterns - identical to AVX2 */
1228 static const uint8_t patterns[8][16] = {
1229 { 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1 },
1230 { 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 },
1231 { 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5 },
1232 { 6, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7 },
1233 { 8, 8, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 9 },
1234 {10,10,10,10,10,10,10,10,11,11,11,11,11,11,11,11 },
1235 {12,12,12,12,12,12,12,12,13,13,13,13,13,13,13,13 },
1236 {14,14,14,14,14,14,14,14,15,15,15,15,15,15,15,15 },
1237 };
1238
1239 const __m128i scale_0 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)patterns[is + 0]));
1240 const __m128i scale_1 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)patterns[is + 1]));
1241 const __m128i scale_2 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)patterns[is + 2]));
1242 const __m128i scale_3 = _mm_shuffle_epi8(scales, _mm_loadu_si128((const __m128i *)patterns[is + 3]));
1243 is += 4;
1244
1245 /* Load low bits */
1246 const __m256i q4bits1 = _mm256_loadu_si256((const __m256i *)q4);
1247 q4 += 32;
1248 const __m256i q4bits2 = _mm256_loadu_si256((const __m256i *)q4);
1249 q4 += 32;
1250 const __m256i q4bitsH = _mm256_loadu_si256((const __m256i *)qh);
1251 qh += 32;
1252
1253 /* Extract high 2-bit contributions */
1254 const __m256i q4h_0 = _mm256_slli_epi16(_mm256_and_si256(q4bitsH, m2), 4);
1255 const __m256i q4h_1 = _mm256_slli_epi16(_mm256_and_si256(_mm256_srli_epi16(q4bitsH, 2), m2), 4);
1256 const __m256i q4h_2 = _mm256_slli_epi16(_mm256_and_si256(_mm256_srli_epi16(q4bitsH, 4), m2), 4);
1257 const __m256i q4h_3 = _mm256_slli_epi16(_mm256_and_si256(_mm256_srli_epi16(q4bitsH, 6), m2), 4);
1258
1259 /* Combine low + high bits to get 6-bit values */
1260 const __m256i q4_0 = _mm256_or_si256(_mm256_and_si256(q4bits1, m4), q4h_0);
1261 const __m256i q4_1 = _mm256_or_si256(_mm256_and_si256(q4bits2, m4), q4h_1);
1262 const __m256i q4_2 = _mm256_or_si256(_mm256_and_si256(_mm256_srli_epi16(q4bits1, 4), m4), q4h_2);
1263 const __m256i q4_3 = _mm256_or_si256(_mm256_and_si256(_mm256_srli_epi16(q4bits2, 4), m4), q4h_3);
1264
1265 /* Load Q8_K values */
1266 const __m256i q8_0 = _mm256_loadu_si256((const __m256i *)q8);
1267 q8 += 32;
1268 const __m256i q8_1 = _mm256_loadu_si256((const __m256i *)q8);
1269 q8 += 32;
1270 const __m256i q8_2 = _mm256_loadu_si256((const __m256i *)q8);
1271 q8 += 32;
1272 const __m256i q8_3 = _mm256_loadu_si256((const __m256i *)q8);
1273 q8 += 32;
1274
1275 /* Compute 32 * q8 for offset */
1276 __m256i q8s_0 = _mm256_maddubs_epi16(m32s, q8_0);
1277 __m256i q8s_1 = _mm256_maddubs_epi16(m32s, q8_1);
1278 __m256i q8s_2 = _mm256_maddubs_epi16(m32s, q8_2);
1279 __m256i q8s_3 = _mm256_maddubs_epi16(m32s, q8_3);
1280
1281 /* Multiply q4 * q8 (unsigned * signed) */
1282 __m256i p16_0 = _mm256_maddubs_epi16(q4_0, q8_0);
1283 __m256i p16_1 = _mm256_maddubs_epi16(q4_1, q8_1);
1284 __m256i p16_2 = _mm256_maddubs_epi16(q4_2, q8_2);
1285 __m256i p16_3 = _mm256_maddubs_epi16(q4_3, q8_3);
1286
1287 /* Subtract offset: (q4 - 32) * q8 = q4*q8 - 32*q8 */
1288 p16_0 = _mm256_sub_epi16(p16_0, q8s_0);
1289 p16_1 = _mm256_sub_epi16(p16_1, q8s_1);
1290 p16_2 = _mm256_sub_epi16(p16_2, q8s_2);
1291 p16_3 = _mm256_sub_epi16(p16_3, q8s_3);
1292
1293 /* Apply scales - produces 8 int32 each */
1294 p16_0 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_0), p16_0);
1295 p16_1 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_1), p16_1);
1296 p16_2 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_2), p16_2);
1297 p16_3 = _mm256_madd_epi16(_mm256_cvtepi8_epi16(scale_3), p16_3);
1298
1299 /* Accumulate all 4 into sumi (same as AVX2) */
1300 sumi = _mm256_add_epi32(sumi, _mm256_add_epi32(p16_0, p16_1));
1301 sumi = _mm256_add_epi32(sumi, _mm256_add_epi32(p16_2, p16_3));
1302 }
1303
1304 /* Scale by d and accumulate */
1305 acc = _mm256_fmadd_ps(_mm256_set1_ps(d), _mm256_cvtepi32_ps(sumi), acc);
1306 }
1307
1308 return ck_q6k_hsum_float_8(acc);
1309}
1310
1311void gemv_q6_k_q8_k_avx512(float *y,
1312 const void *W,
1313 const void *x_q8,
1314 int M, int K)
1315{
1316 if (!y || !W || !x_q8 || M <= 0 || K <= 0) {
1317 return;
1318 }
1319
1320 const block_q6_K *blocks = (const block_q6_K *)W;
1321 const block_q8_K *x = (const block_q8_K *)x_q8;
1322 const int blocks_per_row = K / QK_K;
1323
1324 for (int row = 0; row < M; ++row) {
1325 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
1326 y[row] = dot_q6_k_q8_k_avx512(w_row, x, K);
1327 }
1328}
1329
1330#endif /* __AVX512F__ && __AVX512BW__ */
1331
1332/* ============================================================================
1333 * Dispatch Functions
1334 * ============================================================================ */
1335
1336/**
1337 * @brief Q6_K x Q8_K dot product (single row)
1338 */
1339void vec_dot_q6_k_q8_k(int n, float *s, const void *vx, const void *vy)
1340{
1341 if (!s || !vx || !vy || n <= 0) {
1342 return;
1343 }
1344
1345 const block_q6_K *x = (const block_q6_K *)vx;
1346 const block_q8_K *y = (const block_q8_K *)vy;
1347
1348 /* This is the architecture-neutral scalar oracle. Production x86 dispatch
1349 * uses the separately parity-tested SIMD reduction tree. */
1350 *s = dot_q6_k_q8_k_ref(x, y, n);
1351}
1352
1353/**
1354 * @brief GEMV: y = W @ x where W is Q6_K and x is Q8_K
1355 */
1356void gemv_q6_k_q8_k(float *y,
1357 const void *W,
1358 const void *x_q8,
1359 int M, int K)
1360{
1361 if (ck_strict_parity_enabled() || ck_q6k_q8k_force_ref()) {
1362 gemv_q6_k_q8_k_ref(y, W, x_q8, M, K);
1363 return;
1364 }
1365
1366#if defined(__AVX2__)
1367 /* llama.cpp's x86 Q6_K production graph keeps the AVX2 reduction order
1368 * even when AVX-512 is available. Wider ISA availability is not a license
1369 * to change this numerical contract. */
1370 gemv_q6_k_q8_k_avx2(y, W, x_q8, M, K);
1371 return;
1372#elif defined(__AVX__)
1373 gemv_q6_k_q8_k_avx(y, W, x_q8, M, K);
1374 return;
1375#elif defined(__SSE4_1__)
1376 gemv_q6_k_q8_k_sse(y, W, x_q8, M, K);
1377 return;
1378#elif defined(__ARM_NEON) || defined(__aarch64__)
1379 /* NEON keeps the same integer-exact per-block dot as the scalar
1380 * reference but accumulates the FP32 block contributions sequentially
1381 * instead of through eight partial lanes, so strict parity mode (handled
1382 * above) still uses the reference reduction order. */
1383 gemv_q6_k_q8_k_neon(y, W, x_q8, M, K);
1384 return;
1385#endif
1386 gemv_q6_k_q8_k_ref(y, W, x_q8, M, K);
1387}
1388
1389const char *ck_q6_k_q8_k_provider_name(void)
1390{
1391 if (ck_strict_parity_enabled() || ck_q6k_q8k_force_ref()) {
1392 return "q6_k_q8_k_ref";
1393 }
1394#if defined(__AVX2__)
1395 return "q6_k_q8_k_avx2";
1396#elif defined(__AVX__)
1397 return "q6_k_q8_k_avx";
1398#elif defined(__SSE4_1__)
1399 return "q6_k_q8_k_sse";
1400#elif defined(__ARM_NEON) || defined(__aarch64__)
1401 return "q6_k_q8_k_neon";
1402#else
1403 return "q6_k_q8_k_ref";
1404#endif
1405}
1406
1407/* ============================================================================
1408 * PARALLEL VERSIONS (for parallel orchestration)
1409 *
1410 * These receive ith (thread index) and nth (total threads) from orchestration.
1411 * OpenMP lives in orchestration layer, NOT here.
1412 *
1413 * Naming: *_parallel = receives ith/nth, processes only its portion
1414 * ============================================================================ */
1415
1416/**
1417 * @brief Parallel reference GEMV for Q6_K × Q8_K
1418 *
1419 * Caller provides ith (thread index) and nth (total threads).
1420 * Each thread processes rows [r0, r1).
1421 */
1422void gemv_q6_k_q8_k_parallel(float *y,
1423 const void *W,
1424 const void *x_q8,
1425 int M, int K,
1426 int ith, int nth)
1427{
1428 if (!y || !W || !x_q8 || M <= 0 || K <= 0) return;
1429 if (ith < 0 || nth <= 0 || ith >= nth) return;
1430
1431 /* Compute row range for this thread */
1432 const int dr = (M + nth - 1) / nth;
1433 const int r0 = dr * ith;
1434 const int r1 = (r0 + dr < M) ? (r0 + dr) : M;
1435
1436 if (r0 >= M) return;
1437
1438 const block_q6_K *blocks = (const block_q6_K *)W;
1439 const block_q8_K *x = (const block_q8_K *)x_q8;
1440 const int blocks_per_row = K / QK_K;
1441
1442 for (int row = r0; row < r1; ++row) {
1443 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
1444 y[row] = dot_q6_k_q8_k_ref(w_row, x, K);
1445 }
1446}
1447
1448/**
1449 * @brief Parallel SIMD GEMV for Q6_K × Q8_K
1450 *
1451 * Uses best available SIMD (AVX/SSE) with row prefetching.
1452 * Caller provides ith/nth from OpenMP region.
1453 */
1454void gemv_q6_k_q8_k_parallel_simd(float *y,
1455 const void *W,
1456 const void *x_q8,
1457 int M, int K,
1458 int ith, int nth)
1459{
1460 if (!y || !W || !x_q8 || M <= 0 || K <= 0) return;
1461 if (ith < 0 || nth <= 0 || ith >= nth) return;
1462
1463 const int dr = (M + nth - 1) / nth;
1464 const int r0 = dr * ith;
1465 const int r1 = (r0 + dr < M) ? (r0 + dr) : M;
1466
1467 if (r0 >= M) return;
1468
1469 const block_q6_K *blocks = (const block_q6_K *)W;
1470 const block_q8_K *x = (const block_q8_K *)x_q8;
1471 const int blocks_per_row = K / QK_K;
1472 const int strict = ck_strict_parity_enabled() || ck_q6k_q8k_force_ref();
1473
1474 for (int row = r0; row < r1; ++row) {
1475 const block_q6_K *w_row = blocks + (size_t)row * (size_t)blocks_per_row;
1476#if defined(__AVX2__)
1477 y[row] = strict ? dot_q6_k_q8_k_ref(w_row, x, K)
1478 : dot_q6_k_q8_k_avx2(w_row, x, K);
1479#elif defined(__AVX__)
1480 y[row] = strict ? dot_q6_k_q8_k_ref(w_row, x, K)
1481 : dot_q6_k_q8_k_avx(w_row, x, K);
1482#elif defined(__SSE4_1__)
1483 y[row] = strict ? dot_q6_k_q8_k_ref(w_row, x, K)
1484 : dot_q6_k_q8_k_sse(w_row, x, K);
1485#elif defined(__ARM_NEON) || defined(__aarch64__)
1486 y[row] = strict ? dot_q6_k_q8_k_ref(w_row, x, K)
1487 : dot_q6_k_q8_k_neon(w_row, x, K);
1488#else
1489 y[row] = dot_q6_k_q8_k_ref(w_row, x, K);
1490#endif
1491 }
1492}
1493
1494static inline float ck_dot_q6_k_q8_k_fast_or_ref(const block_q6_K *w,
1495 const block_q8_K *x,
1496 int K);
1497
1498/**
1499 * @brief GEMM: Y = W @ X^T where W is Q6_K and X is Q8_K
1500 *
1501 * @param Y Output matrix [N x M] in row-major
1502 * @param W Weight matrix in Q6_K format [M x K]
1503 * @param X_q8 Input matrix in Q8_K format [N x K]
1504 * @param M Number of output rows (output dim)
1505 * @param N Number of input vectors (batch size)
1506 * @param K Input dimension
1507 */
1508void gemm_q6_k_q8_k(float *Y,
1509 const void *W,
1510 const void *X_q8,
1511 int M, int N, int K)
1512{
1513 if (!Y || !W || !X_q8 || M <= 0 || N <= 0 || K <= 0) {
1514 return;
1515 }
1516
1517 const block_q8_K *X = (const block_q8_K *)X_q8;
1518 const int blocks_per_vec = K / QK_K;
1519
1520 for (int n = 0; n < N; ++n) {
1521 const block_q8_K *x_row = X + (size_t)n * (size_t)blocks_per_vec;
1522 gemv_q6_k_q8_k(&Y[n * M], W, x_row, M, K);
1523 }
1524}
1525
1526/**
1527 * @brief NT GEMM: C = A @ B^T where A is Q8_K and B is Q6_K
1528 *
1529 * This is the typical inference pattern:
1530 * - A: Activations in Q8_K format [M x K]
1531 * - B: Weights in Q6_K format [N x K]
1532 * - C: Output [M x N]
1533 *
1534 * @param A_q8 Input activations in Q8_K format
1535 * @param B Weight matrix in Q6_K format
1536 * @param bias Optional bias vector [N]
1537 * @param C Output matrix
1538 * @param M Batch size (number of tokens)
1539 * @param N Output dimension
1540 * @param K Input dimension
1541 */
1542void gemm_nt_q6_k_q8_k_tiled_parity(const void *A_q8,
1543 const void *B,
1544 const float *bias,
1545 float *C,
1546 int M, int N, int K);
1547
1548void gemm_nt_q6_k_q8_k(const void *A_q8,
1549 const void *B,
1550 const float *bias,
1551 float *C,
1552 int M, int N, int K)
1553{
1554 if (!A_q8 || !B || !C) {
1555 return;
1556 }
1557 if (M <= 0 || N <= 0 || K <= 0) {
1558 return;
1559 }
1560#if defined(__AVX2__) && defined(__FMA__)
1561 if (M >= 8 && ck_strict_parity_enabled()) {
1562 gemm_nt_q6_k_q8_k_tiled_parity(A_q8, B, bias, C, M, N, K);
1563 return;
1564 }
1565#endif
1566
1567 /* Prefill GEMM uses the pinned x86 leaf reduction for short row groups
1568 * and the tiled reduction for strict batches of at least eight rows.
1569 * CK_DEBUG_Q6K_Q8K_REF=1 still forces the scalar reference leaf. */
1570 const block_q8_K *A = (const block_q8_K *)A_q8;
1571 const block_q6_K *W = (const block_q6_K *)B;
1572 const int blocks_per_vec = K / QK_K;
1573 const int blocks_per_row = K / QK_K;
1574
1575 for (int m = 0; m < M; ++m) {
1576 const block_q8_K *a_row = A + (size_t)m * (size_t)blocks_per_vec;
1577 float *c_row = C + (size_t)m * (size_t)N;
1578 for (int n = 0; n < N; ++n) {
1579 const block_q6_K *w_row = W + (size_t)n * (size_t)blocks_per_row;
1580 const float b = bias ? bias[n] : 0.0f;
1581 c_row[n] = ck_dot_q6_k_q8_k_fast_or_ref(w_row, a_row, K) + b;
1582 }
1583 }
1584}
1585
1586/* The pinned x86 GGML tiled graph reduces each 256-value slab in integer
1587 * arithmetic, then applies FP32 scales with one fused accumulation per slab. */
1588void gemm_nt_q6_k_q8_k_tiled_parity_tile(const void *A_q8,
1589 const void *B,
1590 const float *bias,
1591 float *C,
1592 int M, int N, int K,
1593 int m0, int m1, int n0, int n1)
1594{
1595 if (!A_q8 || !B || !C || M <= 0 || N <= 0 || K <= 0 || K % QK_K) return;
1596 if (m0 < 0 || n0 < 0 || m1 > M || n1 > N || m0 >= m1 || n0 >= n1) return;
1597 const block_q8_K *A = (const block_q8_K *) A_q8;
1598 const block_q6_K *W = (const block_q6_K *) B;
1599 const int nb = K / QK_K;
1600 for (int m = m0; m < m1; ++m) {
1601 for (int n = n0; n < n1; ++n) {
1602 float acc = 0.0f;
1603 for (int b = 0; b < nb; ++b) {
1604 const block_q6_K *w = W + (size_t)n * (size_t)nb + b;
1605 const block_q8_K *a = A + (size_t)m * (size_t)nb + b;
1606 int32_t dot = 0;
1607 for (int e = 0; e < QK_K; ++e) {
1608 const int group = (e / 128) * 128;
1609 const int pos = e % 128;
1610 const int lane = pos % 32;
1611 const int ql_off = group / 2 + lane + ((pos / 32) & 1) * 32;
1612 const int qh_off = group / 4 + lane;
1613 const int ql_shift = (pos >= 64) ? 4 : 0;
1614 const int qh_shift = (pos / 32) * 2;
1615 const int code = ((w->ql[ql_off] >> ql_shift) & 15) |
1616 (((w->qh[qh_off] >> qh_shift) & 3) << 4);
1617 dot += (code - 32) * (int)w->scales[e / 16] * (int)a->qs[e];
1618 }
1619 acc = fmaf(GGML_FP16_TO_FP32(w->d) * (float)dot, a->d, acc);
1620 }
1621 C[(size_t)m * (size_t)N + n] = acc + (bias ? bias[n] : 0.0f);
1622 }
1623 }
1624}
1625
1626void gemm_nt_q6_k_q8_k_tiled_parity(const void *A_q8,
1627 const void *B,
1628 const float *bias,
1629 float *C,
1630 int M, int N, int K)
1631{
1632 gemm_nt_q6_k_q8_k_tiled_parity_tile(A_q8, B, bias, C,
1633 M, N, K, 0, M, 0, N);
1634}
1635
1636static void gemm_nt_q6_k_q8_k_prepared_tile_impl(
1637 const void *A_q8,
1638 const void *B_prepared,
1639 const float *bias,
1640 float *C,
1641 int M, int N, int K,
1642 int m0, int m1,
1643 int n0, int n1,
1644 int use_avx512_vnni)
1645{
1646#if !defined(__AVX2__)
1647 (void)A_q8; (void)B_prepared; (void)bias; (void)C;
1648 (void)M; (void)N; (void)K; (void)m0; (void)m1; (void)n0; (void)n1;
1649 (void)use_avx512_vnni;
1650#else
1651#if !defined(__AVX512F__) || !defined(__AVX512BW__) || \
1652 !defined(__AVX512VNNI__)
1653 (void)use_avx512_vnni;
1654#endif
1655 if (!A_q8 || !B_prepared || !C || M <= 0 || N <= 0 || K <= 0 ||
1656 (K % QK_K) != 0) return;
1657 if (m0 < 0) m0 = 0;
1658 if (n0 < 0) n0 = 0;
1659 if (m1 > M) m1 = M;
1660 if (n1 > N) n1 = N;
1661 if (m0 >= m1 || n0 >= n1) return;
1662
1663 const block_q8_K *A = (const block_q8_K *)A_q8;
1664 const block_q6_K_prepared *W =
1665 (const block_q6_K_prepared *)B_prepared;
1666 const int blocks_per_row = K / QK_K;
1667 for (int n = n0; n < n1; ++n) {
1668 const block_q6_K_prepared *w_row =
1669 W + (size_t)n * (size_t)blocks_per_row;
1670 const float b = bias ? bias[n] : 0.0f;
1671 for (int m = m0; m < m1; ++m) {
1672 const block_q8_K *a_row =
1673 A + (size_t)m * (size_t)blocks_per_row;
1674#if defined(__AVX512F__) && defined(__AVX512BW__) && \
1675 defined(__AVX512VNNI__)
1676 const float dot = use_avx512_vnni
1677 ? dot_q6_k_prepared_q8_k_avx512_vnni(w_row, a_row, K)
1678 : dot_q6_k_prepared_q8_k_avx2(w_row, a_row, K);
1679#else
1680 const float dot = dot_q6_k_prepared_q8_k_avx2(w_row, a_row, K);
1681#endif
1682 C[(size_t)m * (size_t)N + (size_t)n] =
1683 dot + b;
1684 }
1685 }
1686#endif
1687}
1688
1689void gemm_nt_q6_k_q8_k_prepared_tile(const void *A_q8,
1690 const void *B_prepared,
1691 const float *bias,
1692 float *C,
1693 int M, int N, int K,
1694 int m0, int m1,
1695 int n0, int n1)
1696{
1697 int use_avx512_vnni = 0;
1698#if defined(__AVX512F__) && defined(__AVX512BW__) && \
1699 defined(__AVX512VNNI__)
1700 use_avx512_vnni = 1;
1701#endif
1702 gemm_nt_q6_k_q8_k_prepared_tile_impl(
1703 A_q8, B_prepared, bias, C, M, N, K,
1704 m0, m1, n0, n1, use_avx512_vnni);
1705}
1706
1707void gemm_nt_q6_k_q8_k_prepared_avx512_vnni(
1708 const void *A_q8,
1709 const void *B_prepared,
1710 const float *bias,
1711 float *C,
1712 int M, int N, int K)
1713{
1714 gemm_nt_q6_k_q8_k_prepared_tile_impl(
1715 A_q8, B_prepared, bias, C, M, N, K,
1716 0, M, 0, N, 1);
1717}
1718
1719void gemm_nt_q6_k_q8_k_prepared(const void *A_q8,
1720 const void *B_prepared,
1721 const float *bias,
1722 float *C,
1723 int M, int N, int K)
1724{
1725 gemm_nt_q6_k_q8_k_prepared_tile(
1726 A_q8, B_prepared, bias, C, M, N, K, 0, M, 0, N);
1727}
1728
1729static inline float ck_dot_q6_k_q8_k_fast_or_ref(const block_q6_K *w,
1730 const block_q8_K *x,
1731 int K)
1732{
1733 if (ck_q6k_q8k_force_ref()) {
1734 return dot_q6_k_q8_k_ref(w, x, K);
1735 }
1736#if defined(__AVX2__)
1737 /* Short strict prefill uses the pinned x86 graph's AVX2 leaf reduction.
1738 * Batched strict prefill uses the separate tiled integer-slab kernel. */
1739 return dot_q6_k_q8_k_avx2(w, x, K);
1740#elif defined(__AVX__)
1741 if (ck_strict_parity_enabled()) return dot_q6_k_q8_k_ref(w, x, K);
1742 return dot_q6_k_q8_k_avx(w, x, K);
1743#elif defined(__SSE4_1__)
1744 if (ck_strict_parity_enabled()) return dot_q6_k_q8_k_ref(w, x, K);
1745 return dot_q6_k_q8_k_sse(w, x, K);
1746#elif defined(__ARM_NEON) || defined(__aarch64__)
1747 if (ck_strict_parity_enabled()) return dot_q6_k_q8_k_ref(w, x, K);
1748 return dot_q6_k_q8_k_neon(w, x, K);
1749#else
1750 return dot_q6_k_q8_k_ref(w, x, K);
1751#endif
1752}
1753
1754/**
1755 * @brief Compute one C[m0:m1, n0:n1] tile for Q8_K activations x Q6_K weights.
1756 *
1757 * Pure tile math only: no threadpool, no global scheduling, no allocation.
1758 * The orchestrator decides how to split tile jobs across cores.
1759 */
1760void gemm_nt_q6_k_q8_k_tile(const void *A_q8,
1761 const void *B,
1762 const float *bias,
1763 float *C,
1764 int M, int N, int K,
1765 int m0, int m1,
1766 int n0, int n1)
1767{
1768 if (!A_q8 || !B || !C) {
1769 return;
1770 }
1771 if (M <= 0 || N <= 0 || K <= 0 || K % QK_K != 0) {
1772 return;
1773 }
1774 if (m0 < 0) m0 = 0;
1775 if (n0 < 0) n0 = 0;
1776 if (m1 > M) m1 = M;
1777 if (n1 > N) n1 = N;
1778 if (m0 >= m1 || n0 >= n1) {
1779 return;
1780 }
1781
1782 const block_q8_K *A = (const block_q8_K *)A_q8;
1783 const block_q6_K *W = (const block_q6_K *)B;
1784 const int blocks_per_vec = K / QK_K;
1785 const int blocks_per_row = K / QK_K;
1786
1787 for (int n = n0; n < n1; ++n) {
1788 const block_q6_K *w_row = W + (size_t)n * (size_t)blocks_per_row;
1789 const float b = bias ? bias[n] : 0.0f;
1790 for (int m = m0; m < m1; ++m) {
1791 const block_q8_K *a_row = A + (size_t)m * (size_t)blocks_per_vec;
1792 C[(size_t)m * (size_t)N + (size_t)n] = ck_dot_q6_k_q8_k_fast_or_ref(w_row, a_row, K) + b;
1793 }
1794 }
1795}
1796
1797void gemm_nt_q6_k_q8_k_m4_tile(const void *A_q8,
1798 const void *B,
1799 const float *bias,
1800 float *C,
1801 int M, int N, int K,
1802 int m0, int m1,
1803 int n0, int n1)
1804{
1805 if (!A_q8 || !B || !C || M <= 0 || N <= 0 || K <= 0 ||
1806 K % QK_K != 0) {
1807 return;
1808 }
1809 if (m0 < 0) m0 = 0;
1810 if (n0 < 0) n0 = 0;
1811 if (m1 > M) m1 = M;
1812 if (n1 > N) n1 = N;
1813 if (m0 >= m1 || n0 >= n1) return;
1814
1815#if defined(__AVX2__)
1816 if (!ck_strict_parity_enabled() && !ck_q6k_q8k_force_ref()) {
1817 const block_q8_K *A = (const block_q8_K *)A_q8;
1818 const block_q6_K *W = (const block_q6_K *)B;
1819 const int blocks_per_vec = K / QK_K;
1820 for (int n = n0; n < n1; ++n) {
1821 const block_q6_K *w_row =
1822 W + (size_t)n * (size_t)blocks_per_vec;
1823 const float b = bias ? bias[n] : 0.0f;
1824 int m = m0;
1825 for (; m + 4 <= m1; m += 4) {
1826 float values[4];
1827 dot_q6_k_q8_k_avx2_m4(
1828 w_row, A + (size_t)m * (size_t)blocks_per_vec,
1829 blocks_per_vec, 4, K, values);
1830 for (int r = 0; r < 4; ++r) {
1831 C[(size_t)(m + r) * (size_t)N + (size_t)n] = values[r] + b;
1832 }
1833 }
1834 if (m < m1) {
1835 float values[4];
1836 const int rows = m1 - m;
1837 dot_q6_k_q8_k_avx2_m4(
1838 w_row, A + (size_t)m * (size_t)blocks_per_vec,
1839 blocks_per_vec, rows, K, values);
1840 for (int r = 0; r < rows; ++r) {
1841 C[(size_t)(m + r) * (size_t)N + (size_t)n] = values[r] + b;
1842 }
1843 }
1844 }
1845 return;
1846 }
1847#endif
1848 gemm_nt_q6_k_q8_k_tile(
1849 A_q8, B, bias, C, M, N, K, m0, m1, n0, n1);
1850}
1851
1852/* Decode M=2 uses the same row-preserving Q6_K/Q8_K arithmetic as prefill.
1853 * The AVX2 tile unpacks each Q6_K weight block once for both sequences. */
1854void gemm_nt_q6_k_q8_k_m2(const void *A_q8, const void *B,
1855 const float *bias, float *C,
1856 int M, int N, int K)
1857{
1858 if (!A_q8 || !B || !C || M != 2 || N <= 0 || K <= 0 ||
1859 (K % QK_K) != 0) return;
1860 gemm_nt_q6_k_q8_k_m4_tile(A_q8, B, bias, C, M, N, K,
1861 0, M, 0, N);
1862}
1863
1864/**
1865 * @brief Experimental single-thread tiled NT GEMM wrapper.
1866 *
1867 * Kept as a separate symbol from gemm_nt_q6_k_q8_k for benchmarks and parity.
1868 * Production prefill should prefer the v8 2D tile scheduler when enabled.
1869 */
1870void gemm_nt_q6_k_q8_k_tiled(const void *A_q8,
1871 const void *B,
1872 const float *bias,
1873 float *C,
1874 int M, int N, int K)
1875{
1876 enum { TILE_M = 8, TILE_N = 16 };
1877 for (int n0 = 0; n0 < N; n0 += TILE_N) {
1878 const int n1 = (n0 + TILE_N < N) ? (n0 + TILE_N) : N;
1879 for (int m0 = 0; m0 < M; m0 += TILE_M) {
1880 const int m1 = (m0 + TILE_M < M) ? (m0 + TILE_M) : M;
1881 gemm_nt_q6_k_q8_k_tile(A_q8, B, bias, C, M, N, K, m0, m1, n0, n1);
1882 }
1883 }
1884}
ckernel_engine.h
ck_strict_parity_enabled
int ck_strict_parity_enabled(void)
Definition ckernel_strict.c:43
ckernel_quant.h
Quantization block structures for weight-only quantization.
GGML_FP16_TO_FP32
#define GGML_FP16_TO_FP32
Definition ckernel_quant.h:451
ggml_half
ck_half ggml_half
Definition ckernel_quant.h:448
QK_K
#define QK_K
Definition ckernel_quant.h:175
gemm_nt_q6_k_q8_k_tile
void gemm_nt_q6_k_q8_k_tile(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K, int m0, int m1, int n0, int n1)
Compute one C[m0:m1, n0:n1] tile for Q8_K activations x Q6_K weights.
Definition gemm_kernels_q6k_q8k.c:1760
gemv_q6_k_q8_k_ref
void gemv_q6_k_q8_k_ref(float *y, const void *W, const void *x_q8, int M, int K)
Definition gemm_kernels_q6k_q8k.c:285
gemv_q6_k_q8_k_sse
void gemv_q6_k_q8_k_sse(float *y, const void *W, const void *x_q8, int M, int K)
ck_q6_k_prepared_block_size
size_t ck_q6_k_prepared_block_size(void)
Definition gemm_kernels_q6k_q8k.c:57
gemv_q6_k_q8_k_avx
void gemv_q6_k_q8_k_avx(float *y, const void *W, const void *x_q8, int M, int K)
ck_q6_k_q8_k_provider_name
const char * ck_q6_k_q8_k_provider_name(void)
Definition gemm_kernels_q6k_q8k.c:1389
gemm_nt_q6_k_q8_k_tiled_parity
void gemm_nt_q6_k_q8_k_tiled_parity(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K)
NT GEMM: C = A @ B^T where A is Q8_K and B is Q6_K.
Definition gemm_kernels_q6k_q8k.c:1626
gemm_nt_q6_k_q8_k_m4_tile
void gemm_nt_q6_k_q8_k_m4_tile(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K, int m0, int m1, int n0, int n1)
Definition gemm_kernels_q6k_q8k.c:1797
ck_dot_q6_k_q8_k_fast_or_ref
static float ck_dot_q6_k_q8_k_fast_or_ref(const block_q6_K *w, const block_q8_K *x, int K)
Definition gemm_kernels_q6k_q8k.c:1729
gemm_q6_k_q8_k
void gemm_q6_k_q8_k(float *Y, const void *W, const void *X_q8, int M, int N, int K)
GEMM: Y = W @ X^T where W is Q6_K and X is Q8_K.
Definition gemm_kernels_q6k_q8k.c:1508
gemm_nt_q6_k_q8_k_prepared_tile_impl
static void gemm_nt_q6_k_q8_k_prepared_tile_impl(const void *A_q8, const void *B_prepared, const float *bias, float *C, int M, int N, int K, int m0, int m1, int n0, int n1, int use_avx512_vnni)
Definition gemm_kernels_q6k_q8k.c:1636
gemv_q6_k_q8_k
void gemv_q6_k_q8_k(float *y, const void *W, const void *x_q8, int M, int K)
GEMV: y = W @ x where W is Q6_K and x is Q8_K.
Definition gemm_kernels_q6k_q8k.c:1356
gemv_q6_k_q8_k_avx2
void gemv_q6_k_q8_k_avx2(float *y, const void *W, const void *x_q8, int M, int K)
vec_dot_q6_k_q8_k
void vec_dot_q6_k_q8_k(int n, float *s, const void *vx, const void *vy)
Q6_K x Q8_K dot product (single row)
Definition gemm_kernels_q6k_q8k.c:1339
gemm_nt_q6_k_q8_k_prepared_avx512_vnni
void gemm_nt_q6_k_q8_k_prepared_avx512_vnni(const void *A_q8, const void *B_prepared, const float *bias, float *C, int M, int N, int K)
Definition gemm_kernels_q6k_q8k.c:1707
gemm_nt_q6_k_q8_k_prepared
void gemm_nt_q6_k_q8_k_prepared(const void *A_q8, const void *B_prepared, const float *bias, float *C, int M, int N, int K)
Definition gemm_kernels_q6k_q8k.c:1719
gemm_nt_q6_k_q8_k_tiled
void gemm_nt_q6_k_q8_k_tiled(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K)
Experimental single-thread tiled NT GEMM wrapper.
Definition gemm_kernels_q6k_q8k.c:1870
gemv_q6_k_q8_k_avx512_vbmi
void gemv_q6_k_q8_k_avx512_vbmi(float *y, const void *W, const void *x_q8, int M, int K)
gemm_nt_q6_k_q8_k_tiled_parity_tile
void gemm_nt_q6_k_q8_k_tiled_parity_tile(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K, int m0, int m1, int n0, int n1)
Definition gemm_kernels_q6k_q8k.c:1588
gemv_q6_k_q8_k_parallel
void gemv_q6_k_q8_k_parallel(float *y, const void *W, const void *x_q8, int M, int K, int ith, int nth)
Parallel reference GEMV for Q6_K × Q8_K.
Definition gemm_kernels_q6k_q8k.c:1422
gemm_nt_q6_k_q8_k_prepared_tile
void gemm_nt_q6_k_q8_k_prepared_tile(const void *A_q8, const void *B_prepared, const float *bias, float *C, int M, int N, int K, int m0, int m1, int n0, int n1)
Definition gemm_kernels_q6k_q8k.c:1689
ck_q6k_q8k_force_ref
static int ck_q6k_q8k_force_ref(void)
Definition gemm_kernels_q6k_q8k.c:108
ck_q6_k_prepare_weight
void ck_q6_k_prepare_weight(const void *src, void *dst, int N, int K)
Definition gemm_kernels_q6k_q8k.c:74
ck_q6_k_prepared_provider_name
const char * ck_q6_k_prepared_provider_name(void)
Definition gemm_kernels_q6k_q8k.c:62
gemm_nt_q6_k_q8_k
void gemm_nt_q6_k_q8_k(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K)
Definition gemm_kernels_q6k_q8k.c:1548
gemm_nt_q6_k_q8_k_m2
void gemm_nt_q6_k_q8_k_m2(const void *A_q8, const void *B, const float *bias, float *C, int M, int N, int K)
Definition gemm_kernels_q6k_q8k.c:1854
gemv_q6_k_q8_k_avx512
void gemv_q6_k_q8_k_avx512(float *y, const void *W, const void *x_q8, int M, int K)
gemv_q6_k_q8_k_parallel_simd
void gemv_q6_k_q8_k_parallel_simd(float *y, const void *W, const void *x_q8, int M, int K, int ith, int nth)
Parallel SIMD GEMV for Q6_K × Q8_K.
Definition gemm_kernels_q6k_q8k.c:1454
dot_q6_k_q8_k_ref
static float dot_q6_k_q8_k_ref(const block_q6_K *w, const block_q8_K *x, int K)
Scalar dot product for Q6_K x Q8_K.
Definition gemm_kernels_q6k_q8k.c:135
C
#define C(color)
Definition show_config.c:39
block_q6_K
Definition ckernel_quant.h:193
block_q6_K::ql
uint8_t ql[256/2]
Definition ckernel_quant.h:194
block_q6_K::scales
int8_t scales[256/16]
Definition ckernel_quant.h:196
block_q6_K::qh
uint8_t qh[256/4]
Definition ckernel_quant.h:195
block_q6_K::d
ck_half d
Definition ckernel_quant.h:197
block_q8_K
Definition ckernel_quant.h:207
block_q8_K::d
float d
Definition ckernel_quant.h:208
block_q8_K::qs
int8_t qs[256]
Definition ckernel_quant.h:209
block_q8_K::bsums
int16_t bsums[256/16]
Definition ckernel_quant.h:210