docs/modules/f16.md

# `bankML/f16.rs` — F16 weights, multiplied as the shipped ggml b11192 multiplies them

## Summary

`f16.rs` (O3, added in 0.3.4) multiplies F16 weight matrices with the bits of the shipped `libggml-cpu-haswell.so`
at b11192. It was read from that tag's `ggml-cpu.c` (`ggml_compute_forward_mul_mat`), `vec.cpp` (`ggml_vec_dot_f16`)
and `llamafile/sgemm.cpp` (`llamafile_sgemm`, `tinyBLAS`).

ggml multiplies an F16 weight two different ways, and the two give different bits for the same row and column:

- **one column** (a decode step, the logits, a one-token micro-batch): `llamafile_sgemm` refuses `n < 2`, so every
  element is `ggml_vec_dot_f16`: four 8-lane accumulators over 32-element steps, a fixed reduction, then any tail
  (`n % 32`) added in double;
- **two columns or more**: `tinyBLAS<8, __m256, …>` when the weight has a multiple of 4 rows and the row a multiple
  of 8; each element is one 8-lane FMA chain over the row, then `hsum`. Otherwise ggml falls back to `vec_dot`.

In both cases the activation is first rounded to f16 (round to nearest even). The module reproduces that choice of
path by shape. The native forward pass uses it for F16 models (SmolLM2-135M-Instruct and mindX's `mindx-gen39`,
both with tied embeddings), and its attention code (`forward.rs`) uses the f16 helpers over the default f16 KV
cache: `vec_dot_u16` for K cache rows, `mad` and `scale` for the V accumulator, `widen` for K and V tiles. The q8_0
KV cache (0.3.9, `BANKML_CACHE_TYPE=q8_0`) does not use them; its kernels are in `forward.rs`
(`attend_head_q8`, `vec_dot_q8_0_q8_0`) and `q1_0.rs` (`quantize_row_q8_0`), see [forward.md](forward.md).

## Technical usage

Constant: `F16_BYTES = 2`.

```rust
pub fn to_f16(x: &[f32]) -> Vec<u16>                     // ggml_cpu_fp32_to_fp16
pub fn dequantize_row(row: &[u8], out: &mut [f32])       // F16 bytes widened to f32 (get_rows)
pub fn widen(src: &[u16], dst: &mut [f32])               // exact; F16C where present
pub fn vec_dot_ref(x: &[u8], y: &[u16]) -> f32           // ggml_vec_dot_f16, the definition
pub fn sgemm_dot_ref(x: &[u8], y: &[u16]) -> f32         // one tinyBLAS element, y.len() % 8 == 0
pub fn vec_dot(x: &[u8], y: &[u16]) -> f32               // vec_dot_ref's bits, fastest path
pub fn vec_dot_u16(x: &[u16], y: &[u16]) -> f32          // the same, weight row as f16 values (a K cache row)
pub fn mad(acc: &mut [u16], v: &[u16], w: f32)           // ggml_vec_mad_f16 (flash attention's V accumulator)
pub fn scale(acc: &mut [u16], s: f32)                    // ggml_vec_scale_f16
pub fn takes_sgemm(rows: usize, k: usize, n: usize) -> bool
```

`x` is an F16 weight row as little-endian bytes; `y` is the f16 activation. `takes_sgemm` is
`n >= 2 && rows % 4 == 0 && k % 8 == 0`.

### Matrix products

```rust
pub fn mat_vec_par(pool: &crate::par::Pool, wb: &[u8], rows: usize, a: &[u16], out: &mut [f32])
pub fn mat_vec_multi_par(pool: &crate::par::Pool, ws: &[(&[u8], usize)], a: &[u16], outs: &mut [&mut [f32]])
pub fn mat_mul_par(pool: &crate::par::Pool, wb: &[u8], rows: usize, cols: &[Vec<u16>], out: &mut [f32])
```

- `mat_vec_par`: one column, `vec_dot` per row, rows in chunks on the pool.
- `mat_vec_multi_par`: several matrices by one column in one pass of the pool (Q, K and V; gate and up). Each output
  row is still one `vec_dot`.
- `mat_mul_par`: writes `out[c·rows + r]`; `sgemm_dot` when `takes_sgemm`, else `vec_dot` per element. On AVX2 the
  sgemm case runs in 4-row by up-to-3-column register blocks (tinyBLAS's own blocking), with the columns widened to
  f32 once per product.

## How it is verified

- `avx2_equals_the_scalar_models`: `vec_dot` and `mat_mul_par` against `vec_dot_ref` / `sgemm_dot_ref` for row
  lengths 8 to 1,536 (tails included), 1–6 columns, 12 and 13 rows (both sides of `takes_sgemm`), every bit.
- `mad_and_scale_match_the_scalar_definitions`: the attention helpers, odd tails included.
- `the_two_reductions_differ`: the two ggml paths give different bits for the same inputs, so the oracle has to
  check the path choice too.
- `oracle_ggml_b11192_f16` (`#[ignore]`, needs `.models/SmolLM2-135M-Instruct-F16.gguf` and `.models/oracle-f16` from
  `testing/f16_oracle.py`): the shipped library computes `ggml_mul_mat` graphs on SmolLM2's real matrices and on
  synthetic shapes that reach the tails and fallbacks. Recorded ([oracles.md](../oracles.md) §5d): 211 of 211 tensors
  widened bit-exact; 552,268 of 552,268 elements in 87 products (61 tinyBLAS, 26 `vec_dot`).
- End to end, the model oracles and llama-server token oracles cover the F16 models (800 of 800 rows each for
  SmolLM2-135M-Instruct and `mindx-gen39`).

## Design notes

Exact float orders, as read from b11192 (changing any of them changes bits):

- **Activation.** F16's `vec_dot_type` is F16 and its `from_float` is `ggml_cpu_fp32_to_fp16` (F16C
  `vcvtps2ph imm=0`, round to nearest even), so the activation is rounded to f16 before any product: `to_f16`.
- **`vec_dot` (`ggml_vec_dot_f16`).** Four 8-lane f32 accumulators fed by FMAs over 32-element steps, reduced
  `(a0 + a2) + (a1 + a3)`, then low half + high half, then two horizontal adds; any tail (`n % 32`) is added in
  double.
- **`sgemm_dot` (`tinyBLAS<8, __m256, __m256, ggml_fp16_t, ggml_fp16_t, float>`).** Taken when the weight has a
  multiple of 4 rows and the row a multiple of 8 (otherwise `llamafile_sgemm` returns false and ggml falls back to
  `vec_dot`). Each element is one 8-lane FMA chain over the row, `hsum`'d as
  `((l0 + l4) + (l2 + l6)) + ((l1 + l5) + (l3 + l7))`.
- **Blocking.** `sgemm_block_avx2` follows tinyBLAS's register blocking: twelve accumulators for 4 rows × 3 columns,
  each step loading four weight vectors and up to three activation vectors (already widened to f32). Every element
  is still its own FMA chain, so the blocking does not change bits.
- **Oracle.** `testing/f16_oracle.py` has the shipped library compute `mul_mat` graphs; every tensor is also widened
  to f32 with ggml's own conversion and compared.

## Advantages and efficiency

- **Bits chosen by shape, as ggml chooses.** Because the two reductions differ, matching ggml requires the same path
  for the same shape. `takes_sgemm` encodes ggml's rule; the tiling and the thread split decide only who computes an
  element, never its bits.
- **F16C and FMA on the hot path.** Weights are widened with `_mm256_cvtph_ps` inside the FMA loop. In `mat_mul_par`
  the activations are widened once per product, not once per weight row.
- **One pool pass for grouped matrices.** `mat_vec_multi_par` hands Q, K and V (or gate and up) to the pool as one
  job, so the pool is woken once instead of per matrix.
- **Measured** (laptop, three threads, [PERFORMANCE.md](../PERFORMANCE.md), 0.3.4): SmolLM2-135M-Instruct F16 decode
  38.0–38.9 tok/s against llama-server's 40.4–42.7 (about 0.93×); an 858-token prompt in 8.6–9.6 s against
  7.8–8.6 s. The steps that got there: the reference f16 dot, scale and accumulate on F16C + FMA (12.8 → 36.1 tok/s
  together with hash-map tensor lookup and norms read once), and tinyBLAS's 4 × 3 blocking for prompts (858-token
  prompt 18.5 → 14.4 s). The rest of 0.3.4's gain was attention work in `forward.rs`. Each change kept the same bits.
- **Rejected with numbers.** A bounded spin before the pool's workers sleep raised `pool.run` from 15.6 µs to
  162.6 µs per call on the 2-core SMT laptop and did not move decode; it was reverted.
- **Practice.** `#[target_feature(enable = "avx2,fma,f16c")]` (or `"avx2,f16c"`) functions with stated safety
  conditions, called only after `q1_0::has_avx2()` and length checks. Reading f16 slices as bytes is guarded by a
  `compile_error!` on big-endian targets.

## Limitations

- Decode on the F16 models is about 0.93× llama-server's (38 vs 41 tok/s on SmolLM2). The pool's condvar hand-off is
  the suspect; a spin was measured and rejected ([TODO.md](../TODO.md)).
- x86_64 AVX2 + FMA + F16C is the only fast path; other CPUs run the scalar models (same bits).
- Little-endian targets only (`compile_error!` otherwise).
- F16 is the first of ggml's standard weight formats done; `Q8_0` and BF16 weight kernels are next, then `Q4_K_M`
  ([TODO.md](../TODO.md), O3). 0.3.9's q8_0 KV cache is not a Q8_0 weight path: a Q8_0 model is still refused.

## See also

- [q1_0.md](q1_0.md) (the f16 ↔ f32 conversions and `has_avx2`), [par.md](par.md), [gguf.md](gguf.md)
- [../oracles.md](../oracles.md) §5d, [../PERFORMANCE.md](../PERFORMANCE.md) (F16 and the Llama graph),
  [../OLLAMA.md](../OLLAMA.md) (O3, O4), [../TODO.md](../TODO.md)