From 3c53ea6810ecfa54459176618f5e04516d136894 Mon Sep 17 00:00:00 2001 From: Michal Harakal Date: Thu, 27 Aug 2026 21:32:14 +0200 Subject: [PATCH 1/2] =?UTF-8?q?feat(#1195):=20thread=20the=20Q4=5FK/Q6=5FK?= =?UTF-8?q?=20matmul=20kernels=20=E2=80=94=20153=20=E2=86=92=2061=20ms/ste?= =?UTF-8?q?p=20on=20the=20Pixel=208a?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Row-partition threading for both orders of both kernels (feed-order and the #1189 row-major variants), through a shared runner (skainet_row_threads). The design is the product of a measured elimination on device, each variant a cooled Qwen2.5-1.5B Q4_K_M decode step: single-threaded (#1190 baseline) 153 ms pthread_create/join per call 994 ms (cpuidle wakeup × ~600/step) parked pool, static quarter chunks 115 ms (big.LITTLE straggle) parked pool, 64-row stealing 227 ms (short scattered windows) parked pool, guided grains 306 ms (still scheduler-starved) spin-then-park pool + guided grains 61 ms ← shipped The decisive piece is the spin: sub-millisecond parallel bursts separated by sleeps never build per-thread utilization, so EAS/schedutil parks workers on little cores at low clocks — the pool was slower than one pegged big core. Workers spin (`yield`) on the job epoch for ~1 ms before parking on the condvar: utilization stays pegged during decode, the scheduler answers with big cores and full clocks, and everyone parks when work stops (no idle burn). Same reason llama.cpp's pool spins. Structure: kernels refactored around shared per-block terms + row-range workers; entries quantize the activation to Q8 once (read-only across threads) and hand the rows to skainet_run_rows. Guided grains off an atomic cursor — long contiguous streams first, small tail last — degrade to an even split on symmetric cores; nothing is tuned to one SoC. Threshold 512 keeps GQA k/v projections (256 rows) single-threaded. MSVC (no pthreads) and any pthread_create failure degrade to the caller's thread. Bit-identity: workers own disjoint out[] ranges and per-row accumulation order never changes. RowMajorMatmulParityTest grows threaded cases with two oracles: threaded-vs-feed-order on permuted bytes, and threaded-vs-1536 independent single-row calls (pins the partition arithmetic itself). Co-authored-by: Claude Fable 5 --- .../native/CMakeLists.txt | 1 + .../native/CMakeLists.txt | 1 + .../native/include/skainet_kernels.h | 14 + .../native/src/q4k_matmul.c | 284 +++++++++--------- .../native/src/q6k_matmul.c | 246 ++++++++------- .../native/src/skainet_row_threads.c | 231 ++++++++++++++ .../native/src/skainet_row_threads.h | 25 ++ .../exec/kernel/RowMajorMatmulParityTest.kt | 41 +++ 8 files changed, 594 insertions(+), 249 deletions(-) create mode 100644 skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.c create mode 100644 skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.h diff --git a/skainet-backends/skainet-backend-jni-cpu/native/CMakeLists.txt b/skainet-backends/skainet-backend-jni-cpu/native/CMakeLists.txt index bcb2be72f..b79dcf28d 100644 --- a/skainet-backends/skainet-backend-jni-cpu/native/CMakeLists.txt +++ b/skainet-backends/skainet-backend-jni-cpu/native/CMakeLists.txt @@ -14,6 +14,7 @@ set(SKAINET_KERNEL_SOURCES ${SKAINET_KERNELS_ROOT}/src/q5_0_matmul.c ${SKAINET_KERNELS_ROOT}/src/q5_1_matmul.c ${SKAINET_KERNELS_ROOT}/src/q8_0_matmul.c + ${SKAINET_KERNELS_ROOT}/src/skainet_row_threads.c ${SKAINET_KERNELS_ROOT}/src/q4k_matmul.c ${SKAINET_KERNELS_ROOT}/src/q5k_matmul.c ${SKAINET_KERNELS_ROOT}/src/q6k_matmul.c diff --git a/skainet-backends/skainet-backend-native-cpu/native/CMakeLists.txt b/skainet-backends/skainet-backend-native-cpu/native/CMakeLists.txt index e97b9d0b0..d6780f0c6 100644 --- a/skainet-backends/skainet-backend-native-cpu/native/CMakeLists.txt +++ b/skainet-backends/skainet-backend-native-cpu/native/CMakeLists.txt @@ -12,6 +12,7 @@ endif() set(SKAINET_KERNEL_SOURCES src/skainet_smoke.c src/skainet_cpu_features.c + src/skainet_row_threads.c src/q4k_matmul.c src/q5k_matmul.c src/q6k_matmul.c diff --git a/skainet-backends/skainet-backend-native-cpu/native/include/skainet_kernels.h b/skainet-backends/skainet-backend-native-cpu/native/include/skainet_kernels.h index 9d4377a15..02ac35694 100644 --- a/skainet-backends/skainet-backend-native-cpu/native/include/skainet_kernels.h +++ b/skainet-backends/skainet-backend-native-cpu/native/include/skainet_kernels.h @@ -48,6 +48,11 @@ SKAINET_API void skainet_smoke_double(const float* input, float* output, int32_t * * Caller owns input/weight/output memory; the kernel does not retain * pointers past return. input_dim must be a multiple of 256. + * + * Threads over output rows (disjoint out[] slices, pthreads, up to 4) + * when output_dim >= 512 (#1195); single-threaded below that and on + * MSVC. Results are bit-identical either way — per output row the + * accumulation order over blocks never changes. */ SKAINET_API void skainet_q4k_matmul( const float* input, @@ -68,6 +73,9 @@ SKAINET_API void skainet_q4k_matmul( * i.e. the bytes exactly as they sit in a .gguf file — which is what lets * mmap'd weights be fed to this kernel with no relayout copy (each row's * blocks are read strictly sequentially). + * + * Threads over output rows when output_dim >= 512 (#1195) — see + * skainet_q4k_matmul; bit-identical to the single-threaded result. */ SKAINET_API void skainet_q4k_matmul_rm( const float* input, @@ -120,6 +128,9 @@ SKAINET_API void skainet_q5k_matmul( * weight + weight_byte_offset + (block_idx * output_dim + o) * 210 * * input_dim must be a multiple of 256. + * + * Threads over output rows when output_dim >= 512 (#1195) — see + * skainet_q4k_matmul; bit-identical to the single-threaded result. */ SKAINET_API void skainet_q6k_matmul( const float* input, @@ -139,6 +150,9 @@ SKAINET_API void skainet_q6k_matmul( * weight + weight_byte_offset + (o * blocks_per_row + block_idx) * 210 * — the bytes exactly as they sit in a .gguf file, so an mmap'd weight * needs no relayout copy. + * + * Threads over output rows when output_dim >= 512 (#1195) — see + * skainet_q4k_matmul; bit-identical to the single-threaded result. */ SKAINET_API void skainet_q6k_matmul_rm( const float* input, diff --git a/skainet-backends/skainet-backend-native-cpu/native/src/q4k_matmul.c b/skainet-backends/skainet-backend-native-cpu/native/src/q4k_matmul.c index fb26d9f75..e5fc011f3 100644 --- a/skainet-backends/skainet-backend-native-cpu/native/src/q4k_matmul.c +++ b/skainet-backends/skainet-backend-native-cpu/native/src/q4k_matmul.c @@ -1,6 +1,7 @@ #include "skainet_kernels.h" #include "skainet_simd.h" #include "skainet_cpu_features.h" +#include "skainet_row_threads.h" #include #include @@ -176,121 +177,119 @@ static void skainet_q4k_block_dot_generic( #endif /* - * Native Q4_K matrix-vector multiply matching the - * sk.ainet.backend.api.kernel.Q4KMatmulKernel SPI contract. Single input row - * times an `outputDim x inputDim` Q4_K-packed weight laid out - * (blockIdx * outputDim + o) * 144 bytes. - * - * Fused int8 dot path (ggml-style): the input row is quantized to Q8 ONCE per - * 256-block (reused across all output rows), then each weight sub-block is an - * int8 dot-product against the Q8 activation: - * acc += d_in[b] * ( d * Σ_s scaleIdx[s]*intDot[s] - dMin * Σ_s minIdx[s]*intSum[s] ) + * One block's contribution to out[o]: + * d_in[b] * ( d * Σ_s scaleIdx[s]*intDot[s] - dMin * Σ_s minIdx[s]*intSum[s] ) * where intDot[s] = Σ q8[i]*code[i] and intSum[s] = Σ q8[i] over the sub-block. * On AArch64 with dotprod (asimddp) the inner dot uses vdotq_s32 (16 int8 MACs * per instruction); otherwise a scalar integer fallback (auto-vectorized). - * The index mapping (groups, lo/hi sub-blocks, input alignment) is identical to - * the previous float kernel, which was parity-checked against Panama. + * Shared by the feed-order and row-major entries — identical math, so the two + * orders (and any row partition, #1195) stay bit-identical per output row. */ -SKAINET_API void skainet_q4k_matmul( - const float* SKAINET_RESTRICT input, - int32_t input_offset, - const uint8_t* SKAINET_RESTRICT weight, - int32_t weight_byte_offset, - int32_t input_dim, - int32_t output_dim, - float* SKAINET_RESTRICT output, - int32_t output_offset +static inline float skainet_q4k_block_term( + const uint8_t* SKAINET_RESTRICT block, + const int8_t* SKAINET_RESTRICT q8_block, + float di, + int use_dp ) { - if (output_dim <= 0 || input_dim <= 0) return; - -#ifdef SKAINET_DOTPROD_DISPATCH - /* One probe per matmul call; cached in skainet_cpu_has_dotprod. */ - const int use_dp = skainet_cpu_has_dotprod(); -#endif - - const int32_t blocks_per_input_dim = input_dim / Q4K_BLOCK_SIZE; - const float* in_base = input + input_offset; - float* out_base = output + output_offset; - - /* Pre-quantize the whole input row to Q8 once (reused across all o). */ - int8_t* q8 = (int8_t*) malloc((size_t) input_dim * sizeof(int8_t)); - float* d_in = (float*) malloc((size_t) blocks_per_input_dim * sizeof(float)); - if (q8 == NULL || d_in == NULL) { free(q8); free(d_in); return; } - for (int32_t b = 0; b < blocks_per_input_dim; ++b) { - d_in[b] = skainet_q8_quantize_block(in_base + (size_t) b * Q4K_BLOCK_SIZE, - q8 + (size_t) b * Q4K_BLOCK_SIZE); - } - int scale_idx[Q4K_SUB_BLOCKS]; int min_idx[Q4K_SUB_BLOCKS]; - /* - * Loop order: block OUTER, output row INNER. The weight is packed - * block-major — (blockIdx * output_dim + o) * 144 — so for a fixed block, - * consecutive `o` are exactly 144 bytes apart: the weight bytes are read - * strictly sequentially (prefetch- and cache-line-friendly). The reverse - * order (o outer) strides output_dim*144 bytes per step (~295 KB on the - * down-proj), which on an in-order A55 with small caches makes every weight - * read a cold miss and dominates runtime regardless of inner-loop compute. - * out_base[o] is accumulated across blocks (output_dim*4 bytes stays hot in - * cache); the accumulation order over blocks is unchanged, so this is - * numerically identical to the o-outer form. - */ - for (int32_t o = 0; o < output_dim; ++o) out_base[o] = 0.0f; + const uint16_t d_bits = (uint16_t) block[0] | ((uint16_t) block[1] << 8); + const uint16_t d_min_bits = (uint16_t) block[2] | ((uint16_t) block[3] << 8); + const float d = skainet_half_to_float(d_bits); + const float d_min = skainet_half_to_float(d_min_bits); - for (int32_t block_idx = 0; block_idx < blocks_per_input_dim; ++block_idx) { - const int8_t* q8_block = q8 + (size_t) block_idx * Q4K_BLOCK_SIZE; - const float di = d_in[block_idx]; - const uint8_t* block = weight + weight_byte_offset - + (size_t)(block_idx * output_dim) * Q4K_BYTES_PER_BLOCK; + skainet_q4k_decode_scales(block + 4, scale_idx, min_idx); - for (int32_t o = 0; o < output_dim; ++o, block += Q4K_BYTES_PER_BLOCK) { - const uint16_t d_bits = (uint16_t) block[0] | ((uint16_t) block[1] << 8); - const uint16_t d_min_bits = (uint16_t) block[2] | ((uint16_t) block[3] << 8); - const float d = skainet_half_to_float(d_bits); - const float d_min = skainet_half_to_float(d_min_bits); + const uint8_t* qs = block + 16; - skainet_q4k_decode_scales(block + 4, scale_idx, min_idx); - - const uint8_t* qs = block + 16; - - int64_t block_scale_dot = 0; - int64_t block_min_sum = 0; + int64_t block_scale_dot = 0; + int64_t block_min_sum = 0; #if defined(SKAINET_HAVE_DOTPROD) - skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); + (void) use_dp; + skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, + &block_scale_dot, &block_min_sum); #elif defined(SKAINET_DOTPROD_DISPATCH) - if (use_dp) { - skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); - } else { - skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); - } + if (use_dp) { + skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, + &block_scale_dot, &block_min_sum); + } else { + skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, + &block_scale_dot, &block_min_sum); + } #else - skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); + (void) use_dp; + skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, + &block_scale_dot, &block_min_sum); #endif - out_base[o] += di * (d * (float) block_scale_dot - d_min * (float) block_min_sum); + return di * (d * (float) block_scale_dot - d_min * (float) block_min_sum); +} + +/* Everything a row-range worker needs; read-only during the parallel section. */ +typedef struct { + const uint8_t* weight_base; /* weight + weight_byte_offset */ + const int8_t* q8; + const float* d_in; + float* out_base; + int32_t blocks_per_input_dim; + int32_t output_dim; + int use_dp; /* meaningful only under SKAINET_DOTPROD_DISPATCH */ +} skainet_q4k_ctx; + +/* + * Feed-order rows [o_start, o_end). Loop order: block OUTER, output row INNER. + * The weight is packed block-major — (blockIdx * output_dim + o) * 144 — so for + * a fixed block this range's rows are one contiguous 144·(o_end−o_start) byte + * run: reads stay sequential (prefetch- and cache-line-friendly; the o-outer + * order would stride output_dim*144 per step, a cold miss per read on an + * in-order A55). out[o] accumulates across blocks in unchanged order, so the + * result is numerically identical to the o-outer form and to any partition. + */ +static void skainet_q4k_rows_feed(void* vctx, int32_t o_start, int32_t o_end) { + const skainet_q4k_ctx* c = (const skainet_q4k_ctx*) vctx; + for (int32_t o = o_start; o < o_end; ++o) c->out_base[o] = 0.0f; + for (int32_t block_idx = 0; block_idx < c->blocks_per_input_dim; ++block_idx) { + const int8_t* q8_block = c->q8 + (size_t) block_idx * Q4K_BLOCK_SIZE; + const float di = c->d_in[block_idx]; + const uint8_t* block = c->weight_base + + ((size_t) block_idx * c->output_dim + o_start) * Q4K_BYTES_PER_BLOCK; + for (int32_t o = o_start; o < o_end; ++o, block += Q4K_BYTES_PER_BLOCK) { + c->out_base[o] += skainet_q4k_block_term(block, q8_block, di, c->use_dp); } } - - free(q8); - free(d_in); } /* - * Row-major variant (#1189): the weight stays in canonical GGUF file order — + * Row-major rows [o_start, o_end) (#1189): canonical GGUF file order — * (o * blocks_per_row + b) * 144 — so an mmap'd tensor is fed as-is, no - * relayout copy. Loop order is o OUTER here: each output row's blocks are - * contiguous on disk, so the weight bytes are still read strictly - * sequentially; the Q8-quantized activation (input_dim bytes) stays hot - * across rows. Per-row accumulation order over blocks matches the feed-order - * kernel's, so results are bit-identical. + * relayout copy. o OUTER: each row's blocks are contiguous on disk, reads stay + * strictly sequential; the Q8 activation (input_dim bytes) stays hot across + * rows. Per-row accumulation order matches the feed-order worker's. */ -SKAINET_API void skainet_q4k_matmul_rm( +static void skainet_q4k_rows_rm(void* vctx, int32_t o_start, int32_t o_end) { + const skainet_q4k_ctx* c = (const skainet_q4k_ctx*) vctx; + const uint8_t* block = c->weight_base + + (size_t) o_start * c->blocks_per_input_dim * Q4K_BYTES_PER_BLOCK; + for (int32_t o = o_start; o < o_end; ++o) { + float acc = 0.0f; + for (int32_t block_idx = 0; block_idx < c->blocks_per_input_dim; + ++block_idx, block += Q4K_BYTES_PER_BLOCK) { + acc += skainet_q4k_block_term(block, c->q8 + (size_t) block_idx * Q4K_BLOCK_SIZE, + c->d_in[block_idx], c->use_dp); + } + c->out_base[o] = acc; + } +} + +/* + * Shared entry: quantize the input row to Q8 once (reused across all rows and + * all threads — read-only after this point), then run the worker over the + * output rows, threaded per skainet_row_threads.h (#1195: ≥512 rows → up to 4 + * pthreads; below that, or on MSVC, the calling thread does all rows). + */ +static void skainet_q4k_matmul_run( const float* SKAINET_RESTRICT input, int32_t input_offset, const uint8_t* SKAINET_RESTRICT weight, @@ -298,19 +297,14 @@ SKAINET_API void skainet_q4k_matmul_rm( int32_t input_dim, int32_t output_dim, float* SKAINET_RESTRICT output, - int32_t output_offset + int32_t output_offset, + skainet_row_range_fn worker ) { if (output_dim <= 0 || input_dim <= 0) return; -#ifdef SKAINET_DOTPROD_DISPATCH - const int use_dp = skainet_cpu_has_dotprod(); -#endif - const int32_t blocks_per_input_dim = input_dim / Q4K_BLOCK_SIZE; const float* in_base = input + input_offset; - float* out_base = output + output_offset; - /* Pre-quantize the whole input row to Q8 once (reused across all o). */ int8_t* q8 = (int8_t*) malloc((size_t) input_dim * sizeof(int8_t)); float* d_in = (float*) malloc((size_t) blocks_per_input_dim * sizeof(float)); if (q8 == NULL || d_in == NULL) { free(q8); free(d_in); return; } @@ -319,50 +313,64 @@ SKAINET_API void skainet_q4k_matmul_rm( q8 + (size_t) b * Q4K_BLOCK_SIZE); } - int scale_idx[Q4K_SUB_BLOCKS]; - int min_idx[Q4K_SUB_BLOCKS]; - - const uint8_t* block = weight + weight_byte_offset; - for (int32_t o = 0; o < output_dim; ++o) { - float acc = 0.0f; - for (int32_t block_idx = 0; block_idx < blocks_per_input_dim; - ++block_idx, block += Q4K_BYTES_PER_BLOCK) { - const int8_t* q8_block = q8 + (size_t) block_idx * Q4K_BLOCK_SIZE; - const float di = d_in[block_idx]; - - const uint16_t d_bits = (uint16_t) block[0] | ((uint16_t) block[1] << 8); - const uint16_t d_min_bits = (uint16_t) block[2] | ((uint16_t) block[3] << 8); - const float d = skainet_half_to_float(d_bits); - const float d_min = skainet_half_to_float(d_min_bits); - - skainet_q4k_decode_scales(block + 4, scale_idx, min_idx); - - const uint8_t* qs = block + 16; - - int64_t block_scale_dot = 0; - int64_t block_min_sum = 0; - -#if defined(SKAINET_HAVE_DOTPROD) - skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); -#elif defined(SKAINET_DOTPROD_DISPATCH) - if (use_dp) { - skainet_q4k_block_dot_dp(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); - } else { - skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); - } + skainet_q4k_ctx ctx; + ctx.weight_base = weight + weight_byte_offset; + ctx.q8 = q8; + ctx.d_in = d_in; + ctx.out_base = output + output_offset; + ctx.blocks_per_input_dim = blocks_per_input_dim; + ctx.output_dim = output_dim; +#ifdef SKAINET_DOTPROD_DISPATCH + /* One probe per matmul call; cached in skainet_cpu_has_dotprod. */ + ctx.use_dp = skainet_cpu_has_dotprod(); #else - skainet_q4k_block_dot_generic(qs, q8_block, scale_idx, min_idx, - &block_scale_dot, &block_min_sum); + ctx.use_dp = 0; #endif - acc += di * (d * (float) block_scale_dot - d_min * (float) block_min_sum); - } - out_base[o] = acc; - } + skainet_run_rows(worker, &ctx, output_dim); free(q8); free(d_in); } + +/* + * Native Q4_K matrix-vector multiply matching the + * sk.ainet.backend.api.kernel.Q4KMatmulKernel SPI contract. Single input row + * times an `outputDim x inputDim` Q4_K-packed weight laid out + * (blockIdx * outputDim + o) * 144 bytes. Threads over output rows when + * outputDim >= 512 (#1195); parity-checked against Panama. + */ +SKAINET_API void skainet_q4k_matmul( + const float* SKAINET_RESTRICT input, + int32_t input_offset, + const uint8_t* SKAINET_RESTRICT weight, + int32_t weight_byte_offset, + int32_t input_dim, + int32_t output_dim, + float* SKAINET_RESTRICT output, + int32_t output_offset +) { + skainet_q4k_matmul_run(input, input_offset, weight, weight_byte_offset, + input_dim, output_dim, output, output_offset, + skainet_q4k_rows_feed); +} + +/* + * Row-major variant (#1189): the weight stays in canonical GGUF file order — + * (o * blocks_per_row + b) * 144 — see skainet_q4k_rows_rm. Bit-identical to + * the feed-order kernel; threads over output rows when outputDim >= 512. + */ +SKAINET_API void skainet_q4k_matmul_rm( + const float* SKAINET_RESTRICT input, + int32_t input_offset, + const uint8_t* SKAINET_RESTRICT weight, + int32_t weight_byte_offset, + int32_t input_dim, + int32_t output_dim, + float* SKAINET_RESTRICT output, + int32_t output_offset +) { + skainet_q4k_matmul_run(input, input_offset, weight, weight_byte_offset, + input_dim, output_dim, output, output_offset, + skainet_q4k_rows_rm); +} diff --git a/skainet-backends/skainet-backend-native-cpu/native/src/q6k_matmul.c b/skainet-backends/skainet-backend-native-cpu/native/src/q6k_matmul.c index cb8d46150..ea7d78742 100644 --- a/skainet-backends/skainet-backend-native-cpu/native/src/q6k_matmul.c +++ b/skainet-backends/skainet-backend-native-cpu/native/src/q6k_matmul.c @@ -1,6 +1,7 @@ #include "skainet_kernels.h" #include "skainet_simd.h" #include "skainet_cpu_features.h" +#include "skainet_row_threads.h" #include #include @@ -152,18 +153,100 @@ static int64_t skainet_q6k_weighted_dot_generic(const int8_t* SKAINET_RESTRICT q #endif /* - * Native Q6_K matrix-vector multiply matching the - * sk.ainet.backend.api.kernel.Q6KMatmulKernel SPI contract. A single - * input row times an `outputDim x inputDim` Q6_K-packed weight tensor - * laid out (blockIdx * outputDim + o) * 210 bytes. - * - * Fused int8 dot path (ggml-style, mirrors q4k_matmul.c): the input row is - * quantized to Q8 ONCE per 256-block (reused across all output rows), the 6-bit - * weight is unpacked to centered int8 codes, and each scale-group is an int8 - * dot (vdotq_s32 on dotprod targets) — no 256-float scratch, no per-element - * float multiply. acc = d · d_in · Σ_g sc[g]·Σ_{i∈g} q8[i]·codes[i]. + * One block's contribution to out[o]: the 6-bit weight is unpacked to centered + * int8 codes and each scale-group is an int8 dot (vdotq_s32 on dotprod + * targets) — acc term = d · d_in · Σ_g sc[g]·Σ_{i∈g} q8[i]·codes[i]. + * `codes` is the caller's per-thread 256-byte scratch (#1195). Shared by the + * feed-order and row-major entries, so both orders (and any row partition) + * stay bit-identical per output row. */ -SKAINET_API void skainet_q6k_matmul( +static inline float skainet_q6k_block_term( + const uint8_t* SKAINET_RESTRICT block, + const int8_t* SKAINET_RESTRICT q8_block, + float di, + int use_dp, + int8_t* SKAINET_RESTRICT codes +) { + const uint16_t d_bits = (uint16_t) block[Q6K_D_OFFSET] + | ((uint16_t) block[Q6K_D_OFFSET + 1] << 8); + const float d = skainet_q6k_half_to_float(d_bits); + const int8_t* sc = (const int8_t*)(block + Q6K_SCALES_OFFSET); + + skainet_q6k_unpack_codes(block, codes); +#if defined(SKAINET_HAVE_DOTPROD) + (void) use_dp; + const int64_t wdot = skainet_q6k_weighted_dot_dp(q8_block, codes, sc); +#elif defined(SKAINET_DOTPROD_DISPATCH) + const int64_t wdot = use_dp + ? skainet_q6k_weighted_dot_dp(q8_block, codes, sc) + : skainet_q6k_weighted_dot_generic(q8_block, codes, sc); +#else + (void) use_dp; + const int64_t wdot = skainet_q6k_weighted_dot_generic(q8_block, codes, sc); +#endif + + return d * di * (float) wdot; +} + +/* Everything a row-range worker needs; read-only during the parallel section. */ +typedef struct { + const uint8_t* weight_base; /* weight + weight_byte_offset */ + const int8_t* q8; + const float* d_in; + float* out_base; + int32_t blocks_per_input_dim; + int32_t output_dim; + int use_dp; /* meaningful only under SKAINET_DOTPROD_DISPATCH */ +} skainet_q6k_ctx; + +/* + * Feed-order rows [o_start, o_end): block OUTER, row INNER — see q4k_matmul.c + * for the cache rationale. This range's rows are one contiguous run per block; + * out[o] accumulates across blocks in unchanged order under any partition. + */ +static void skainet_q6k_rows_feed(void* vctx, int32_t o_start, int32_t o_end) { + const skainet_q6k_ctx* c = (const skainet_q6k_ctx*) vctx; + int8_t codes[Q6K_BLOCK_SIZE]; + for (int32_t o = o_start; o < o_end; ++o) c->out_base[o] = 0.0f; + for (int32_t block_idx = 0; block_idx < c->blocks_per_input_dim; ++block_idx) { + const int8_t* q8_block = c->q8 + (size_t) block_idx * Q6K_BLOCK_SIZE; + const float di = c->d_in[block_idx]; + const uint8_t* block = c->weight_base + + ((size_t) block_idx * c->output_dim + o_start) * Q6K_BYTES_PER_BLOCK; + for (int32_t o = o_start; o < o_end; ++o, block += Q6K_BYTES_PER_BLOCK) { + c->out_base[o] += skainet_q6k_block_term(block, q8_block, di, c->use_dp, codes); + } + } +} + +/* + * Row-major rows [o_start, o_end) (#1189): canonical GGUF file order — + * (o * blocks_per_row + b) * 210 — an mmap'd tensor is fed as-is, no relayout + * copy; each row's blocks are contiguous on disk. Per-row accumulation order + * matches the feed-order worker's. + */ +static void skainet_q6k_rows_rm(void* vctx, int32_t o_start, int32_t o_end) { + const skainet_q6k_ctx* c = (const skainet_q6k_ctx*) vctx; + int8_t codes[Q6K_BLOCK_SIZE]; + const uint8_t* block = c->weight_base + + (size_t) o_start * c->blocks_per_input_dim * Q6K_BYTES_PER_BLOCK; + for (int32_t o = o_start; o < o_end; ++o) { + float acc = 0.0f; + for (int32_t block_idx = 0; block_idx < c->blocks_per_input_dim; + ++block_idx, block += Q6K_BYTES_PER_BLOCK) { + acc += skainet_q6k_block_term(block, c->q8 + (size_t) block_idx * Q6K_BLOCK_SIZE, + c->d_in[block_idx], c->use_dp, codes); + } + c->out_base[o] = acc; + } +} + +/* + * Shared entry: quantize the input row to Q8 once (read-only afterwards, shared + * by all threads), then run the worker over the output rows, threaded per + * skainet_row_threads.h (#1195). + */ +static void skainet_q6k_matmul_run( const float* SKAINET_RESTRICT input, int32_t input_offset, const uint8_t* SKAINET_RESTRICT weight, @@ -171,20 +254,14 @@ SKAINET_API void skainet_q6k_matmul( int32_t input_dim, int32_t output_dim, float* SKAINET_RESTRICT output, - int32_t output_offset + int32_t output_offset, + skainet_row_range_fn worker ) { if (output_dim <= 0 || input_dim <= 0) return; -#ifdef SKAINET_DOTPROD_DISPATCH - /* One probe per matmul call; cached in skainet_cpu_has_dotprod. */ - const int use_dp = skainet_cpu_has_dotprod(); -#endif - const int32_t blocks_per_input_dim = input_dim / Q6K_BLOCK_SIZE; const float* in_base = input + input_offset; - float* out_base = output + output_offset; - /* Pre-quantize the whole input row to Q8 once (reused across all o). */ int8_t* q8 = (int8_t*) malloc((size_t) input_dim * sizeof(int8_t)); float* d_in = (float*) malloc((size_t) blocks_per_input_dim * sizeof(float)); if (q8 == NULL || d_in == NULL) { free(q8); free(d_in); return; } @@ -193,56 +270,52 @@ SKAINET_API void skainet_q6k_matmul( q8 + (size_t) b * Q6K_BLOCK_SIZE); } - int8_t codes[Q6K_BLOCK_SIZE]; - - /* - * Loop order: block OUTER, output row INNER — see q4k_matmul.c for the - * rationale. The weight is block-major (blockIdx*output_dim + o)*210, so for - * a fixed block consecutive `o` are 210 bytes apart: the weight bytes are - * read sequentially (cache/prefetch friendly) instead of striding - * output_dim*210 per step. out_base[o] accumulates across blocks; the order - * over blocks is unchanged. - */ - for (int32_t o = 0; o < output_dim; ++o) out_base[o] = 0.0f; - - for (int32_t block_idx = 0; block_idx < blocks_per_input_dim; ++block_idx) { - const int8_t* q8_block = q8 + (size_t) block_idx * Q6K_BLOCK_SIZE; - const float di = d_in[block_idx]; - const uint8_t* block = weight + weight_byte_offset - + (size_t)(block_idx * output_dim) * Q6K_BYTES_PER_BLOCK; - - for (int32_t o = 0; o < output_dim; ++o, block += Q6K_BYTES_PER_BLOCK) { - const uint16_t d_bits = (uint16_t) block[Q6K_D_OFFSET] - | ((uint16_t) block[Q6K_D_OFFSET + 1] << 8); - const float d = skainet_q6k_half_to_float(d_bits); - const int8_t* sc = (const int8_t*)(block + Q6K_SCALES_OFFSET); - - skainet_q6k_unpack_codes(block, codes); -#if defined(SKAINET_HAVE_DOTPROD) - const int64_t wdot = skainet_q6k_weighted_dot_dp(q8_block, codes, sc); -#elif defined(SKAINET_DOTPROD_DISPATCH) - const int64_t wdot = use_dp - ? skainet_q6k_weighted_dot_dp(q8_block, codes, sc) - : skainet_q6k_weighted_dot_generic(q8_block, codes, sc); + skainet_q6k_ctx ctx; + ctx.weight_base = weight + weight_byte_offset; + ctx.q8 = q8; + ctx.d_in = d_in; + ctx.out_base = output + output_offset; + ctx.blocks_per_input_dim = blocks_per_input_dim; + ctx.output_dim = output_dim; +#ifdef SKAINET_DOTPROD_DISPATCH + /* One probe per matmul call; cached in skainet_cpu_has_dotprod. */ + ctx.use_dp = skainet_cpu_has_dotprod(); #else - const int64_t wdot = skainet_q6k_weighted_dot_generic(q8_block, codes, sc); + ctx.use_dp = 0; #endif - out_base[o] += d * di * (float) wdot; - } - } + skainet_run_rows(worker, &ctx, output_dim); free(q8); free(d_in); } +/* + * Native Q6_K matrix-vector multiply matching the + * sk.ainet.backend.api.kernel.Q6KMatmulKernel SPI contract. A single + * input row times an `outputDim x inputDim` Q6_K-packed weight tensor + * laid out (blockIdx * outputDim + o) * 210 bytes. Threads over output + * rows when outputDim >= 512 (#1195). + */ +SKAINET_API void skainet_q6k_matmul( + const float* SKAINET_RESTRICT input, + int32_t input_offset, + const uint8_t* SKAINET_RESTRICT weight, + int32_t weight_byte_offset, + int32_t input_dim, + int32_t output_dim, + float* SKAINET_RESTRICT output, + int32_t output_offset +) { + skainet_q6k_matmul_run(input, input_offset, weight, weight_byte_offset, + input_dim, output_dim, output, output_offset, + skainet_q6k_rows_feed); +} + /* * Row-major variant (#1189): the weight stays in canonical GGUF file order — - * (o * blocks_per_row + b) * 210 — so an mmap'd tensor is fed as-is, no - * relayout copy. o OUTER: each row's blocks are contiguous on disk, so the - * weight bytes are still read sequentially; the Q8 activation stays hot. - * Per-row accumulation order over blocks matches the feed-order kernel's, - * so results are bit-identical. + * (o * blocks_per_row + b) * 210 — see skainet_q6k_rows_rm. Bit-identical to + * the feed-order kernel; threads over output rows when outputDim >= 512. */ SKAINET_API void skainet_q6k_matmul_rm( const float* SKAINET_RESTRICT input, @@ -254,56 +327,7 @@ SKAINET_API void skainet_q6k_matmul_rm( float* SKAINET_RESTRICT output, int32_t output_offset ) { - if (output_dim <= 0 || input_dim <= 0) return; - -#ifdef SKAINET_DOTPROD_DISPATCH - const int use_dp = skainet_cpu_has_dotprod(); -#endif - - const int32_t blocks_per_input_dim = input_dim / Q6K_BLOCK_SIZE; - const float* in_base = input + input_offset; - float* out_base = output + output_offset; - - /* Pre-quantize the whole input row to Q8 once (reused across all o). */ - int8_t* q8 = (int8_t*) malloc((size_t) input_dim * sizeof(int8_t)); - float* d_in = (float*) malloc((size_t) blocks_per_input_dim * sizeof(float)); - if (q8 == NULL || d_in == NULL) { free(q8); free(d_in); return; } - for (int32_t b = 0; b < blocks_per_input_dim; ++b) { - d_in[b] = skainet_q6k_q8_quantize_block(in_base + (size_t) b * Q6K_BLOCK_SIZE, - q8 + (size_t) b * Q6K_BLOCK_SIZE); - } - - int8_t codes[Q6K_BLOCK_SIZE]; - - const uint8_t* block = weight + weight_byte_offset; - for (int32_t o = 0; o < output_dim; ++o) { - float acc = 0.0f; - for (int32_t block_idx = 0; block_idx < blocks_per_input_dim; - ++block_idx, block += Q6K_BYTES_PER_BLOCK) { - const int8_t* q8_block = q8 + (size_t) block_idx * Q6K_BLOCK_SIZE; - const float di = d_in[block_idx]; - - const uint16_t d_bits = (uint16_t) block[Q6K_D_OFFSET] - | ((uint16_t) block[Q6K_D_OFFSET + 1] << 8); - const float d = skainet_q6k_half_to_float(d_bits); - const int8_t* sc = (const int8_t*)(block + Q6K_SCALES_OFFSET); - - skainet_q6k_unpack_codes(block, codes); -#if defined(SKAINET_HAVE_DOTPROD) - const int64_t wdot = skainet_q6k_weighted_dot_dp(q8_block, codes, sc); -#elif defined(SKAINET_DOTPROD_DISPATCH) - const int64_t wdot = use_dp - ? skainet_q6k_weighted_dot_dp(q8_block, codes, sc) - : skainet_q6k_weighted_dot_generic(q8_block, codes, sc); -#else - const int64_t wdot = skainet_q6k_weighted_dot_generic(q8_block, codes, sc); -#endif - - acc += d * di * (float) wdot; - } - out_base[o] = acc; - } - - free(q8); - free(d_in); + skainet_q6k_matmul_run(input, input_offset, weight, weight_byte_offset, + input_dim, output_dim, output, output_offset, + skainet_q6k_rows_rm); } diff --git a/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.c b/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.c new file mode 100644 index 000000000..10a319ca8 --- /dev/null +++ b/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.c @@ -0,0 +1,231 @@ +/* + * Row-range threading shared by the packed matmul kernels (#1195). + * + * The packed kernels parallelize over OUTPUT ROWS: participants pull disjoint + * `[o_start, o_end)` grains of the output vector, so there are no + * accumulation races and no synchronization beyond the completion barrier. + * Per output row the accumulation order over input blocks is unchanged, and + * a row's result does not depend on which thread computes it — threaded + * results are bit-identical to single-threaded ones (the parity suites are + * the oracle for that claim). + * + * Execution model, shaped by a series of Pixel 8a measurements (#1195, all + * Qwen2.5-1.5B Q4_K_M decode steps; single-threaded baseline 153 ms): + * + * 1. POOL, NOT CREATE/JOIN. The obvious create/join-per-call variant ran + * 994 ms/step: ~600 pthread_create per step against cores in deep cpuidle + * pays milliseconds of wakeup latency each. Workers are created once, + * lazily (pthread_once), and live for the process. + * + * 2. SPIN BRIEFLY, THEN PARK. A pool whose workers sleep between jobs still + * measured only 115–306 ms/step across chunking variants — sub-millisecond + * parallel bursts separated by sleeps never accumulate per-thread + * utilization, so EAS/schedutil keeps the workers on little cores at low + * clocks while the single-threaded caller would have pegged one big core + * at max clock. Workers therefore spin (`yield`) on the job epoch for + * ~SKAINET_SPIN_ITERS before parking on the condvar: during a decode the + * gaps between matmuls are far shorter than the spin window, utilization + * stays pegged, and the scheduler answers with big cores and full clocks + * (the same reason llama.cpp's thread pool spins). Once work stops + * arriving, everyone parks — no battery burn at idle. + * + * 3. GUIDED grains — remaining/(2·parts), floored at SKAINET_MATMUL_GRAIN — + * off an atomic cursor: long contiguous streams first (prefetch-friendly), + * shrinking toward the tail so a straggler holds a small tail rather than + * a quarter of the matrix. On symmetric cores this degrades to an even + * split; nothing here is tuned to one SoC's topology. + * + * Threading engages only when `n` reaches the threshold — below it fixed + * costs dominate any win, and tiny projections (e.g. GQA k/v with output_dim + * 256) stay single-threaded on purpose. A concurrent second caller while the + * pool is busy simply runs its rows on its own thread (correct, unshared); + * so does everything if worker creation ever failed. + * + * MSVC has no ; there the runner degrades to a plain call on the + * caller's thread. + */ +#include "skainet_row_threads.h" + +#if !defined(_MSC_VER) + +#include +#include + +#define SKAINET_POOL_WORKERS (SKAINET_MATMUL_THREADS - 1) +#define SKAINET_MATMUL_GRAIN 64 +/* ~a millisecond of `yield`s — longer than the gaps between a decode step's + * matmul calls, far shorter than "the model stopped decoding". */ +#define SKAINET_SPIN_ITERS (1u << 20) + +#if defined(__aarch64__) || defined(__arm__) +#define SKAINET_CPU_RELAX() __asm__ __volatile__("yield" ::: "memory") +#elif defined(__x86_64__) || defined(__i386__) +#define SKAINET_CPU_RELAX() __asm__ __volatile__("pause" ::: "memory") +#else +#define SKAINET_CPU_RELAX() ((void) 0) +#endif + +typedef struct { + pthread_mutex_t m; /* serializes callers; guards the park/wake handoffs */ + pthread_cond_t cv_work; + pthread_cond_t cv_done; + skainet_row_range_fn fn; /* job fields: written before the epoch release-store */ + void* ctx; + int32_t n; + int parts; + atomic_int_fast32_t cursor; /* next unclaimed row of the current job */ + atomic_ulong epoch; /* release-published per job; the workers' work signal */ + atomic_int remaining; /* workers not yet finished with the current job */ + atomic_int work_waiters; /* workers parked on cv_work (broadcast only then) */ + atomic_int done_waiter; /* caller parked on cv_done (signal only then) */ + int busy; /* a job is in flight (second callers go solo) */ + int workers_alive; +} skainet_row_pool; + +static skainet_row_pool skainet_g_row_pool = { + PTHREAD_MUTEX_INITIALIZER, PTHREAD_COND_INITIALIZER, PTHREAD_COND_INITIALIZER, + NULL, NULL, 0, 0, 0, 0UL, 0, 0, 0, 0, 0, +}; +static pthread_once_t skainet_g_row_pool_once = PTHREAD_ONCE_INIT; + +/* Guided sizing (OpenMP `schedule(guided)` shape): see file header, point 3. */ +static void skainet_row_pool_drain(skainet_row_range_fn fn, void* ctx, int32_t n, int parts) { + int_fast32_t cur = atomic_load_explicit(&skainet_g_row_pool.cursor, memory_order_relaxed); + for (;;) { + if ((int32_t) cur >= n) return; + int32_t want = (n - (int32_t) cur) / (2 * parts); + if (want < SKAINET_MATMUL_GRAIN) want = SKAINET_MATMUL_GRAIN; + if (atomic_compare_exchange_weak_explicit( + &skainet_g_row_pool.cursor, &cur, cur + want, + memory_order_relaxed, memory_order_relaxed)) { + const int32_t s = (int32_t) cur; + int32_t e = s + want; + if (e > n) e = n; + fn(ctx, s, e); + cur = atomic_load_explicit(&skainet_g_row_pool.cursor, memory_order_relaxed); + } + /* CAS failure reloaded `cur`; loop retries with the fresh value. */ + } +} + +static void* skainet_row_pool_worker(void* arg) { + (void) arg; + unsigned long seen = 0UL; + for (;;) { + /* Spin for the next job; park only when none arrives in the window. */ + unsigned spins = 0; + while (atomic_load_explicit(&skainet_g_row_pool.epoch, memory_order_acquire) == seen) { + if (++spins >= SKAINET_SPIN_ITERS) { + pthread_mutex_lock(&skainet_g_row_pool.m); + atomic_fetch_add_explicit(&skainet_g_row_pool.work_waiters, 1, memory_order_relaxed); + while (atomic_load_explicit(&skainet_g_row_pool.epoch, memory_order_acquire) == seen) { + pthread_cond_wait(&skainet_g_row_pool.cv_work, &skainet_g_row_pool.m); + } + atomic_fetch_sub_explicit(&skainet_g_row_pool.work_waiters, 1, memory_order_relaxed); + pthread_mutex_unlock(&skainet_g_row_pool.m); + break; + } + SKAINET_CPU_RELAX(); + } + seen = atomic_load_explicit(&skainet_g_row_pool.epoch, memory_order_acquire); + + /* Job fields were written before the epoch release-store — the acquire + * above orders these plain reads, and they are stable while busy. */ + skainet_row_pool_drain(skainet_g_row_pool.fn, skainet_g_row_pool.ctx, + skainet_g_row_pool.n, skainet_g_row_pool.parts); + + if (atomic_fetch_sub_explicit(&skainet_g_row_pool.remaining, 1, memory_order_acq_rel) == 1) { + /* Last one out: the park/wake handoff must be decided under the + * mutex — a bare done_waiter load could miss a caller that is + * between setting the flag and blocking, and sleep it forever. + * Cost: one lock/unlock per job, by one thread. */ + pthread_mutex_lock(&skainet_g_row_pool.m); + if (atomic_load_explicit(&skainet_g_row_pool.done_waiter, memory_order_relaxed) != 0) { + pthread_cond_signal(&skainet_g_row_pool.cv_done); + } + pthread_mutex_unlock(&skainet_g_row_pool.m); + } + } + /* unreachable */ +} + +static void skainet_row_pool_init(void) { + for (int i = 0; i < SKAINET_POOL_WORKERS; ++i) { + pthread_t t; + if (pthread_create(&t, NULL, skainet_row_pool_worker, NULL) != 0) { + break; /* fewer participants; the cursor still covers every row */ + } + pthread_detach(t); + ++skainet_g_row_pool.workers_alive; + } +} + +/* + * Run `fn` over rows [0, n): every participant (workers + the calling + * thread) pulls guided grains from the cursor; the caller then spins briefly + * on the completion count before parking. Single-threaded when n is under + * the threshold, when the pool is busy with another caller's job, or when no + * workers exist. + */ +void skainet_run_rows(skainet_row_range_fn fn, void* ctx, int32_t n) { + if (n <= 0) return; + if (n < SKAINET_MATMUL_THREAD_THRESHOLD) { + fn(ctx, 0, n); + return; + } + pthread_once(&skainet_g_row_pool_once, skainet_row_pool_init); + + pthread_mutex_lock(&skainet_g_row_pool.m); + if (skainet_g_row_pool.workers_alive == 0 || skainet_g_row_pool.busy) { + pthread_mutex_unlock(&skainet_g_row_pool.m); + fn(ctx, 0, n); + return; + } + skainet_g_row_pool.busy = 1; + skainet_g_row_pool.fn = fn; + skainet_g_row_pool.ctx = ctx; + skainet_g_row_pool.n = n; + skainet_g_row_pool.parts = skainet_g_row_pool.workers_alive + 1; + atomic_store_explicit(&skainet_g_row_pool.cursor, 0, memory_order_relaxed); + atomic_store_explicit(&skainet_g_row_pool.remaining, skainet_g_row_pool.workers_alive, + memory_order_relaxed); + /* Publish: job fields above happen-before this release-store. */ + atomic_store_explicit(&skainet_g_row_pool.epoch, + atomic_load_explicit(&skainet_g_row_pool.epoch, memory_order_relaxed) + 1UL, + memory_order_release); + if (atomic_load_explicit(&skainet_g_row_pool.work_waiters, memory_order_relaxed) != 0) { + pthread_cond_broadcast(&skainet_g_row_pool.cv_work); + } + const int parts = skainet_g_row_pool.parts; + pthread_mutex_unlock(&skainet_g_row_pool.m); + + skainet_row_pool_drain(fn, ctx, n, parts); + + /* Spin briefly for the stragglers' tail, then park. */ + unsigned spins = 0; + while (atomic_load_explicit(&skainet_g_row_pool.remaining, memory_order_acquire) > 0) { + if (++spins >= SKAINET_SPIN_ITERS) { + pthread_mutex_lock(&skainet_g_row_pool.m); + atomic_store_explicit(&skainet_g_row_pool.done_waiter, 1, memory_order_release); + while (atomic_load_explicit(&skainet_g_row_pool.remaining, memory_order_acquire) > 0) { + pthread_cond_wait(&skainet_g_row_pool.cv_done, &skainet_g_row_pool.m); + } + atomic_store_explicit(&skainet_g_row_pool.done_waiter, 0, memory_order_relaxed); + pthread_mutex_unlock(&skainet_g_row_pool.m); + break; + } + SKAINET_CPU_RELAX(); + } + + pthread_mutex_lock(&skainet_g_row_pool.m); + skainet_g_row_pool.busy = 0; + pthread_mutex_unlock(&skainet_g_row_pool.m); +} + +#else /* MSVC: no pthreads — single-threaded, same numerics */ + +void skainet_run_rows(skainet_row_range_fn fn, void* ctx, int32_t n) { + if (n > 0) fn(ctx, 0, n); +} + +#endif diff --git a/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.h b/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.h new file mode 100644 index 000000000..6f78f57d6 --- /dev/null +++ b/skainet-backends/skainet-backend-native-cpu/native/src/skainet_row_threads.h @@ -0,0 +1,25 @@ +/* + * Row-range threading shared by the packed matmul kernels (#1195) — see + * skainet_row_threads.c for the pool and the measured rationale. + */ +#ifndef SKAINET_ROW_THREADS_H +#define SKAINET_ROW_THREADS_H + +#include + +#define SKAINET_MATMUL_THREADS 4 +#define SKAINET_MATMUL_THREAD_THRESHOLD 512 + +typedef void (*skainet_row_range_fn)(void* ctx, int32_t o_start, int32_t o_end); + +/* + * Run `fn` over rows [0, n): threaded over the shared worker pool when + * n >= SKAINET_MATMUL_THREAD_THRESHOLD, on the calling thread otherwise + * (also on MSVC, when the pool is busy with another caller, or when worker + * creation failed). Bit-identical to the single-threaded result: workers own + * disjoint [o_start, o_end) row ranges and per-row accumulation order never + * changes. + */ +void skainet_run_rows(skainet_row_range_fn fn, void* ctx, int32_t n); + +#endif /* SKAINET_ROW_THREADS_H */ diff --git a/skainet-backends/skainet-backend-native-cpu/src/jvmTest/kotlin/sk/ainet/exec/kernel/RowMajorMatmulParityTest.kt b/skainet-backends/skainet-backend-native-cpu/src/jvmTest/kotlin/sk/ainet/exec/kernel/RowMajorMatmulParityTest.kt index e9fd850f2..5ec8da6d6 100644 --- a/skainet-backends/skainet-backend-native-cpu/src/jvmTest/kotlin/sk/ainet/exec/kernel/RowMajorMatmulParityTest.kt +++ b/skainet-backends/skainet-backend-native-cpu/src/jvmTest/kotlin/sk/ainet/exec/kernel/RowMajorMatmulParityTest.kt @@ -148,4 +148,45 @@ class RowMajorMatmulParityTest { @Test fun q6k_single_block_multi_row() = assertQ6kParity(256, 16, seed = 8) @Test fun q6k_multi_block_multi_row() = assertQ6kParity(1024, 64, seed = 321) @Test fun q6k_honors_weight_byte_offset() = assertQ6kParity(512, 8, seed = 18, pad = 129) + + // ---- #1195: outputDim >= 512 engages the row-partition threading. Two oracles: ---- + + /** Threaded full-matrix call vs the feed-order kernel on permuted bytes (both threaded). */ + @Test fun q4k_threaded_parity_vs_feed_order() = assertQ4kParity(512, 1536, seed = 77) + @Test fun q6k_threaded_parity_vs_feed_order() = assertQ6kParity(512, 1536, seed = 78) + + /** + * Threaded full-matrix call vs 1536 independent single-row calls (each under the + * threshold, so single-threaded) — pins the partition arithmetic itself: every row of a + * threaded call must be bit-identical to that row computed alone. + */ + @Test + fun q4k_threaded_equals_per_row_calls() { + val inputDim = 512 + val n = 1536 + val bpr = inputDim / BLOCK + val rowMajor = randomBlocks(bpr * n, Q4K_BPB, intArrayOf(0, 2), seed = 91) + val input = FloatArray(inputDim) { Random(91 + it).nextFloat() - 0.5f } + + val full = callRm(q4kRm!!, input, rowMajor, 0, inputDim, n) + for (o in 0 until n step 97) { + val single = callRm(q4kRm!!, input, rowMajor, o * bpr * Q4K_BPB, inputDim, 1) + assertEquals(single[0].toRawBits(), full[o].toRawBits(), "Q4_K row $o: threaded diverged from solo") + } + } + + @Test + fun q6k_threaded_equals_per_row_calls() { + val inputDim = 512 + val n = 1536 + val bpr = inputDim / BLOCK + val rowMajor = randomBlocks(bpr * n, Q6K_BPB, intArrayOf(208), seed = 92) + val input = FloatArray(inputDim) { Random(92 + it).nextFloat() - 0.5f } + + val full = callRm(q6kRm!!, input, rowMajor, 0, inputDim, n) + for (o in 0 until n step 97) { + val single = callRm(q6kRm!!, input, rowMajor, o * bpr * Q6K_BPB, inputDim, 1) + assertEquals(single[0].toRawBits(), full[o].toRawBits(), "Q6_K row $o: threaded diverged from solo") + } + } } From f281e2fe95b4ea00b7f1f6329c99a5d1418421ac Mon Sep 17 00:00:00 2001 From: Michal Harakal Date: Thu, 27 Aug 2026 21:59:22 +0200 Subject: [PATCH 2/2] harness(#1193): residency knob for the M2-A5 harness + the mixed-quant fallback trap, measured MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `-e residency heap` restores heap staging (default mapped) so staging strategies can be A/B'd on models that fit the cap. First use found a trap instead of an answer: SmolLM2-135M's hidden size (576) is not a 256-multiple, so llama.cpp quantized most of its matrices as Q8_0 (k-quant fallback) — and heap Q8_0 in canonical order without prepack has no BLOCKED_ROW_MAJOR kernel, so dispatch silently served the decoding reference: heap + prepack=true 66 ms/step mapped + prepack=false 48,771 ms/step (Q8_0 -> reference, ~800x) mapped + prepack=true 65 ms/step (mapped rm + prepacked Q8_0) Exactly #1193's "silent fallback" case and #1192's Q8_0 priority, now with numbers. The Qwen runs never hit it (all dims 256-multiples, pure Q4_K/Q6_K). Co-authored-by: Claude Fable 5 --- .../ainet/exec/harness/android/M2A5DeviceMeasurement.kt | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/skainet-backends/skainet-backend-jni-cpu/src/androidTest/kotlin/sk/ainet/exec/harness/android/M2A5DeviceMeasurement.kt b/skainet-backends/skainet-backend-jni-cpu/src/androidTest/kotlin/sk/ainet/exec/harness/android/M2A5DeviceMeasurement.kt index a1d382494..44145bf76 100644 --- a/skainet-backends/skainet-backend-jni-cpu/src/androidTest/kotlin/sk/ainet/exec/harness/android/M2A5DeviceMeasurement.kt +++ b/skainet-backends/skainet-backend-jni-cpu/src/androidTest/kotlin/sk/ainet/exec/harness/android/M2A5DeviceMeasurement.kt @@ -105,6 +105,10 @@ class M2A5DeviceMeasurement { val prepack = args.getString("prepack")?.toBoolean() ?: false val steps = arg("steps", 16) val warmup = arg("warmup", 4) + // `residency=heap` (default mapped) restores heap staging — the #1193 A/B lever for + // models that fit the cap. The plan below prices the same form the load uses (#1190). + val residency = if (args.getString("residency") == "heap") WeightResidency.HEAP else WeightResidency.MAPPED + val loadForm = WeightForm(shape = WeightShapeOrientation.OUT_IN, residency = residency) val report = StringBuilder() fun line(s: String = "") { report.append(s).append('\n') } @@ -116,13 +120,10 @@ class M2A5DeviceMeasurement { line("- device: ${Build.MANUFACTURER} ${Build.MODEL}, Android ${Build.VERSION.RELEASE} (SDK ${Build.VERSION.SDK_INT}), ABI ${Build.SUPPORTED_ABIS.firstOrNull()}") line("- ART heap cap (Runtime.maxMemory): ${mb(heapCap)}") line("- model: ${modelFile.name}, ${mb(modelFile.length())} on disk") - line("- ctx=$ctxLen, decode steps=$steps (warm-up $warmup), prepack=$prepack") + line("- ctx=$ctxLen, decode steps=$steps (warm-up $warmup), prepack=$prepack, residency=${residency.name.lowercase()}") line() // ---- plan, from the header only -------------------------------------------------- - // The plan gets the same WeightForm the load below uses, so mapped-servable weights are - // budgeted against the page cache instead of the heap cap (#1189). - val loadForm = WeightForm(shape = WeightShapeOrientation.OUT_IN, residency = WeightResidency.MAPPED) val (plan, geometry) = MappedRandomAccessSource.open(modelPath).let { src -> StreamingGGUFReader.open(src).use { reader -> val input = reader.planInput(ctx = ctxLen, formFor = { loadForm })