Skip to content

et-backend: F16 vecdot GEMV + matrix-engine GEMM #28

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
RehanQasim-dev wants to merge 5 commits into
base: et
Choose a base branch
from
+1,341 −257
Draft

et-backend: F16 vecdot GEMV + matrix-engine GEMM #28

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
31265d3
ggml-et: centralize scp_signal/scp_wait in platform.h
RehanQasim-dev Jul 30, 2026
2fdd04b
et-backend: F16 GEMV rewrite for register-resident vecdot
RehanQasim-dev Jul 14, 2026
f3690db
et-backend: F16 matrix-engine GEMM kernel for prefill
RehanQasim-dev Jul 23, 2026
0f5fcbf
ET mul_mat: add F16-weight/F32-activation matrix-engine kernel
RehanQasim-dev Jul 30, 2026
f1bea73
ggml-et: clarify F16 vecdot fallback comment covers two thresholds
RehanQasim-dev Jul 30, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions ggml/src/ggml-et/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ set(KERNELS
mul_mat_Q4_0_matrix_engine
mul_mat_f16
mul_mat_f16_matrix_engine
mul_mat_f16_f32_matrix_engine
rope_f32
unary_f32
sqr_f32
Expand Down
105 changes: 105 additions & 0 deletions ggml/src/ggml-et/et-kernels/src/block_ops.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -670,6 +670,111 @@ static inline float compute_block_dot_product_f16(const uint16_t * a_block, cons
return compute_block_dot_product_f16_partial(a_block, b_col_start, QK_F16);
}

//******************************************************************************
// Hoisted F16 row-dot API (register-resident accumulator)
//
// Mirrors the q8_dot_* API above: the running sum lives in vector registers
// for the whole row instead of being spilled to memory every 8 elements, so
// the inner loop is load + convert + fmadd with no accumulator round-trip.
//
// Register contract (shared convention with q8_dot_*):
// f20-f23 row accumulators (persistent across tiles, reset per row)
// f31 gather pattern (byte offsets of 8 consecutive f16)
// f11-f18 scratch within tile
// t0, f1-f5 scratch within reduce
// Caller sets the vector mask to 0xFF once around the surrounding loop.
//******************************************************************************

static inline void __attribute__((always_inline))
f16_dot_reset(void) {
// Four independent lane accumulators break the fmadd dependency chain so
// multiple A/B loads stay in flight, hiding load latency.
__asm__ volatile(
"fbci.pi f20, 0\n"
"fbci.pi f21, 0\n"
"fbci.pi f22, 0\n"
"fbci.pi f23, 0\n"
::: "f20", "f21", "f22", "f23");
}

// Accumulate n_blocks blocks of QK_F16 (=32) f16 A values times f32 B into
// f20..f23 (one accumulator per 8-lane chunk within the block).
static inline void __attribute__((always_inline))
f16_dot_tile(const uint16_t * a_row, const float * b_col, int64_t n_blocks) {
static const int32_t gather_pattern[8] = {0, 2, 4, 6, 8, 10, 12, 14};
__asm__ volatile("flw.ps f31, %[g]\n"
: : [g] "m"(*(const int32_t(*)[8])gather_pattern)
: "f31");

// fgh.ps/flw.ps are blocking loads, so a non-blocking prefetch hint issued
// a few blocks ahead turns a DRAM miss into an L1 hit by the time the
// gather runs.
const int64_t PF = 4;
for (int64_t kb = 0; kb < n_blocks; kb++) {
const uint16_t * a_ptr = a_row + (kb << 5); // 32 f16 per block
const float * b_ptr = b_col + (kb << 5);
if (kb + PF < n_blocks) {
const uint16_t * pfa = a_ptr + (PF << 5); // A block = 64B = 1 line
const float * pfb = b_ptr + (PF << 5); // B block = 128B = 2 lines
__asm__ volatile(
"lb x0, 0(%[pa])\n"
"lb x0, 0(%[pb])\n"
"lb x0, 64(%[pb])\n"
: : [pa] "r"(pfa), [pb] "r"(pfb));
}
__asm__ volatile(
// Issue all 8 independent memory ops first (4 A gathers + 4 B loads)
// so several misses are outstanding before any consumer stalls.
"fgh.ps f11, f31(%[a0])\n"
"fgh.ps f13, f31(%[a1])\n"
"fgh.ps f15, f31(%[a2])\n"
"fgh.ps f17, f31(%[a3])\n"
"flw.ps f12, %[b0]\n"
"flw.ps f14, %[b1]\n"
"flw.ps f16, %[b2]\n"
"flw.ps f18, %[b3]\n"
"fcvt.ps.f16 f11, f11\n"
"fcvt.ps.f16 f13, f13\n"
"fcvt.ps.f16 f15, f15\n"
"fcvt.ps.f16 f17, f17\n"
"fmadd.ps f20, f11, f12, f20\n"
"fmadd.ps f21, f13, f14, f21\n"
"fmadd.ps f22, f15, f16, f22\n"
"fmadd.ps f23, f17, f18, f23\n"
:
: [a0] "r"(a_ptr), [a1] "r"(a_ptr + 8),
[a2] "r"(a_ptr + 16), [a3] "r"(a_ptr + 24),
[b0] "m"(*(const float(*)[8])&b_ptr[0]),
[b1] "m"(*(const float(*)[8])&b_ptr[8]),
[b2] "m"(*(const float(*)[8])&b_ptr[16]),
[b3] "m"(*(const float(*)[8])&b_ptr[24])
: "f11", "f12", "f13", "f14", "f15", "f16", "f17", "f18",
"f20", "f21", "f22", "f23"
);
}
}

static inline float __attribute__((always_inline))
f16_dot_reduce(void) {
float result;
__asm__ __volatile__ (
// Combine the 4 lane accumulators, then horizontal-sum the 8 lanes.
"fadd.ps f20, f20, f21, rne \n\t"
"fadd.ps f22, f22, f23, rne \n\t"
"fadd.ps f20, f20, f22, rne \n\t"
"fswizz.ps f1, f20, 0xB1 \n\t"
"fadd.ps f2, f20, f1, rne \n\t"
"fswizz.ps f3, f2, 0x4E \n\t"
"fadd.ps f4, f2, f3, rne \n\t"
"fmvz.x.ps t0, f4, 4 \n\t"
"fbcx.ps f5, t0 \n\t"
"fadd.ps %[vout], f4, f5, rne \n\t"
: [vout] "=f" (result)
:: "t0", "f1", "f2", "f3", "f4", "f5", "f20", "f21", "f22", "f23"
);
return result;
}

// Compute dot product between f32 block and f32 column vector
// Vectorized: processes 8 elements at a time using ET vector instructions
// Block size: up to 16 f32 values (can handle partial blocks for misaligned K)
Expand Down
19 changes: 0 additions & 19 deletions ggml/src/ggml-et/et-kernels/src/mul_mat_Q4_0_matrix_engine.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -33,25 +33,6 @@
#define SCP_CONSUMED_OFF (SCP_READY_OFF + 64) // 4160
#define SCP_PER_MINION (SCP_CONSUMED_OFF + 64) // 4224

// Signal a counter value to the other hart via L2 SCP.
static inline void __attribute__((always_inline)) scp_signal(volatile uint32_t * flag, uint32_t value) {
*flag = value;
FENCE;
evict_to_l2((const void *) flag, 1, 64);
WAIT_CACHEOPS;
}

// Wait for a counter in L2 SCP to reach the expected value.
static inline void __attribute__((always_inline)) scp_wait(volatile uint32_t * flag, uint32_t expected) {
while (1) {
evict_to_l2((const void *) flag, 1, 64);
WAIT_CACHEOPS;
if (*flag >= expected) {
return;
}
}
}

// Dequantize one 32-element Q4_0 block of TILE_M weight rows into the FP32
// panel, written directly in TenB [k][m] order: panel[k*TILE_M + m].
// Low nibble of byte i -> k = i
Expand Down
175 changes: 96 additions & 79 deletions ggml/src/ggml-et/et-kernels/src/mul_mat_f16.c
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,13 +1,21 @@
//******************************************************************************
// MUL_MAT Kernel
// MUL_MAT Kernel (F16 weights)
// Matrix multiplication: C[M,N] = A[M,K] * B[K,N]
//
// Decode-optimized distribution: one output element (row m of A dotted with
// column n of B) is the unit of work, striped by thread_id across every hart
// of every active shire. This keeps all 32 shires busy for a GEMV
// (M=4096, N=1) instead of leaving most idle, and the register-resident
// f16_dot row helper keeps the accumulator in registers for the whole K
// dimension.
//******************************************************************************

#include "block_ops.h"
#include "ggml_tensor.h"
#include "math_fp.h"
#include "platform.h"
#include "quants.h"
#include "tensor.h"

#include <stdint.h>

Expand All @@ -18,18 +26,13 @@ int entry_point(struct ggml_et_binary_params * params, void * env) {
return -1;
}

// Thread coordination
// Thread coordination, use every hart of every active shire.
int thread_id = get_relative_thread_id(kernel_env->shire_mask);
int num_threads = get_num_threads(kernel_env->shire_mask);

if (thread_id < 0 || (thread_id & 1)) {
return 0; // Skip odd threads to avoid resource contention
if (thread_id < 0) {
return 0;
}

int effective_thread_id = thread_id / 2;
int effective_num_threads = (num_threads + 1) / 2;

// Extract tensor references
struct ggml_tensor * src0 = &params->src0; // Weight matrix A (F16)
struct ggml_tensor * src1 = &params->src1; // Activation matrix B (F16/F32)
struct ggml_tensor * dst = &params->dst; // Output matrix C (F32)
Expand All @@ -43,7 +46,7 @@ int entry_point(struct ggml_et_binary_params * params, void * env) {
const uint16_t * src0_data = (const uint16_t *) src0->data;
float * dst_data = (float *) dst->data;

// Dimensions and Strides
// Dimensions and strides
const int64_t K = src0->ne[0];
const int64_t M = src0->ne[1];
const int64_t N = src1->ne[1];
Expand All @@ -56,87 +59,101 @@ int entry_point(struct ggml_et_binary_params * params, void * env) {
const size_t nb11 = src1->nb[1], nb12 = src1->nb[2], nb13 = src1->nb[3];
const size_t nb1 = dst->nb[1], nb2 = dst->nb[2], nb3 = dst->nb[3];

// F16 specific block size (Usually QK_F16)
const int block_size = QK_F16;
const int block_size = QK_F16; // 32
const int64_t K_blocks = K / block_size;
const int64_t K_remainder = K % block_size;

// Threading distribution
const uint64_t total_elements = M * N * ne2 * ne3;
const uint64_t per_thread = 16;
const uint64_t threads_stride = per_thread * effective_num_threads;

if (effective_thread_id * per_thread >= total_elements) {
return 0;
}

// Broadcasting support
const int64_t r2 = ne12 / ne02;
const int64_t r3 = ne13 / ne03;

for (uint64_t base_idx = effective_thread_id * per_thread; base_idx < total_elements; base_idx += threads_stride) {
for (uint64_t j = 0; j < per_thread; j++) {
const uint64_t idx = base_idx + j;
if (idx >= total_elements) {
break;
const int is_f32_b = (src1->type == GGML_TYPE_F32);
const uint64_t total_elements = (uint64_t) M * N * ne2 * ne3;

// Stage the reused B activation vector into per-shire L2 SCP. Streaming the
// weight matrix thrashes L2 and evicts B, so B is re-read from DRAM
// repeatedly during the decode GEMV. Stage it once per shire into L2 SCP (a
// separate SRAM partition cache streaming cannot evict) via
// et_tensor_load_l2scp; one hart per shire issues the DMA loop, then a
// shire barrier lets all harts read B on-chip. Only applies to the common
// decode case: F32 B, a single contiguous non-broadcast B that fits the SCP
// budget.
const int b_contig = (nb11 == (size_t) K * sizeof(float));
const uint64_t b_lines = ((uint64_t) N * K * sizeof(float) + 63) / 64;
const int stage_b = is_f32_b && b_contig &&
ne12 == 1 && ne13 == 1 && ne02 == 1 && ne03 == 1 &&
b_lines <= 8192; // <= 512 KB, within the SCP budget
const float * b_scp = (const float *) et_shire_l2scp_local(0);

if (stage_b) {
if ((get_hart_id() & 63) == 0) {
et_tensor_load_l2scp_conf_t conf;
conf.use_tmask = false;
conf.stride = 64; // advance one 64B cache line per line loaded
uint64_t remaining = b_lines;
uint64_t dst = 0;
uint64_t addr = (uint64_t) src1->data;
while (remaining > 0) {
uint64_t cl = (remaining >= 16) ? 16 : remaining;
conf.dst_start = dst;
conf.addr = addr;
conf.num_lines = cl - 1; // 4-bit field encodes (lines - 1)
conf.id = 0;
et_tensor_load_l2scp(&conf);
WAIT_TENSOR_LOAD_L2_0;
dst += cl;
addr += cl * 64;
remaining -= cl;
}
}
et_barrier(ET_BARRIER_SHIRE); // B is now resident in L2 SCP for all harts
}

// Index decoding
const int64_t i3 = idx / (M * N * ne2);
const int64_t rem3 = idx % (M * N * ne2);
const int64_t i2 = rem3 / (M * N);
const int64_t rem2 = rem3 % (M * N);
const int64_t n = rem2 / M;
const int64_t m = rem2 % M;

const int64_t i03 = i3 / r3, i02 = i2 / r2;
const int64_t i13 = (ne13 > 1) ? i3 : 0, i12 = (ne12 > 1) ? i2 : 0;

float sum = 0.0f;
const uint16_t * f16_row =
(const uint16_t *) ((const char *) src0_data + m * nb01 + i02 * nb02 + i03 * nb03);

if (src1->type == GGML_TYPE_F32) {
const float * src1_data = (const float *) src1->data;

for (int64_t kb = 0; kb < K_blocks; kb++) {
const float * b_col_ptr =
(const float *) ((const char *) src1_data + (kb * block_size) * sizeof(float) + n * nb11 +
i12 * nb12 + i13 * nb13);
sum += compute_block_dot_product_f16_naive(&f16_row[kb * block_size], b_col_ptr);
}

if (K_remainder > 0) {
const int64_t offset = K_blocks * block_size;
const float * b_col_ptr = (const float *) ((const char *) src1_data + offset * sizeof(float) +
n * nb11 + i12 * nb12 + i13 * nb13);
sum += compute_block_dot_product_f16_partial(&f16_row[offset], b_col_ptr, K_remainder);
}
} else {
const uint16_t * src1_data = (const uint16_t *) src1->data;

for (int64_t kb = 0; kb < K_blocks; kb++) {
const uint16_t * b_col_ptr =
(const uint16_t *) ((const char *) src1_data + (kb * block_size) * sizeof(uint16_t) + n * nb11 +
i12 * nb12 + i13 * nb13);
sum += compute_block_dot_product_f16_f16_partial(&f16_row[kb * block_size], b_col_ptr, block_size);
}

if (K_remainder > 0) {
const int64_t offset = K_blocks * block_size;
const uint16_t * b_col_ptr =
(const uint16_t *) ((const char *) src1_data + offset * sizeof(uint16_t) + n * nb11 +
i12 * nb12 + i13 * nb13);
sum += compute_block_dot_product_f16_f16_partial(&f16_row[offset], b_col_ptr, K_remainder);
}
// Set the vector mask (all 8 lanes) once for the whole row loop.
unsigned long saved_mask;
__asm__ volatile("mova.x.m %0" : "=r"(saved_mask));
__asm__ volatile("mov.m.x m0, x0, 0xFF");

for (uint64_t idx = (uint64_t) thread_id; idx < total_elements; idx += (uint64_t) num_threads) {
// Index decoding
const int64_t i3 = idx / (M * N * ne2);
const int64_t rem3 = idx % (M * N * ne2);
const int64_t i2 = rem3 / (M * N);
const int64_t rem2 = rem3 % (M * N);
const int64_t n = rem2 / M;
const int64_t m = rem2 % M;

const int64_t i03 = i3 / r3, i02 = i2 / r2;
const int64_t i13 = (ne13 > 1) ? i3 : 0, i12 = (ne12 > 1) ? i2 : 0;

const uint16_t * f16_row =
(const uint16_t *) ((const char *) src0_data + m * nb01 + i02 * nb02 + i03 * nb03);

float sum;
if (is_f32_b) {
const float * b_col = stage_b
? (b_scp + n * K)
: (const float *) ((const char *) src1->data + n * nb11 + i12 * nb12 + i13 * nb13);

f16_dot_reset();
f16_dot_tile(f16_row, b_col, K_blocks);
sum = f16_dot_reduce();

if (K_remainder > 0) {
const int64_t offset = K_blocks * block_size;
sum += compute_block_dot_product_f16_partial(&f16_row[offset], &b_col[offset], K_remainder);
}

// Atomic store for output
volatile float * c_element =
(volatile float *) ((char *) dst_data + m * dst->nb[0] + n * nb1 + i2 * nb2 + i3 * nb3);
atomic_store_f32(c_element, sum);
} else {
const uint16_t * b_col =
(const uint16_t *) ((const char *) src1->data + n * nb11 + i12 * nb12 + i13 * nb13);
sum = compute_block_dot_product_f16_f16_partial(f16_row, b_col, K);
}

volatile float * c_element =
(volatile float *) ((char *) dst_data + m * dst->nb[0] + n * nb1 + i2 * nb2 + i3 * nb3);
atomic_store_f32(c_element, sum);
}

__asm__ volatile("mova.m.x %0" ::"r"(saved_mask));
return 0;
}
Loading