TL;DR — Building a blockwise INT4 inference engine from scratch forces you to confront the full stack: quantization theory, memory layout optimization, SIMD-friendly kernel design, and fused operator implementation. This guide walks you through every step — packed-nibble weight storage, per-group scale/zero-point computation, and a fused dequantizing GEMM kernel — with real, runnable code you can ship to a GitHub repo and point to in interviews.
Quantization-aware inference engines are the backbone of every major LLM serving stack today. From Meta’s LLM.int8() to NVIDIA’s TensorRT-LLM and the llama.cpp ecosystem, efficient low-precision inference is what makes large models run on consumer hardware. Building even a minimal version of this stack yourself is one of the most signal-rich projects an engineer can add to a portfolio — it demonstrates simultaneous fluency in numerical methods, memory architecture, and high-performance code.
This guide gives you a complete, step-by-step path to build a blockwise INT4 inference engine. You will implement packed-nibble weight storage, per-group quantization parameters, and — the centerpiece — a fused dequantizing GEMM kernel that avoids ever materializing full-precision weights in memory. Every code snippet is real and runnable.
Why This Project Stands Out on a CV
Hiring managers and senior engineers scan portfolios for projects that demonstrate depth across layers, not just surface-level familiarity. This project signals several things simultaneously:
- Systems-level thinking: You’re managing memory layouts (nibble packing), cache behavior, and vectorized instructions — the kind of details that separate someone who has used PyTorch from someone who understands what happens under the hood.
- Numerical computing literacy: Understanding how scale and zero-point interact with rounding, overflow, and accumulation is a skill directly relevant to compiler optimization, embedded systems, and accelerator design.
- Performance engineering: Fusing two operations (dequantize + GEMM) into a single kernel is a canonical optimization pattern that appears in every high-performance stack, from database query execution to graphics pipelines.
- ML infrastructure depth: This sits at the intersection of ML and systems — exactly the intersection where the most interesting engineering jobs live (infra ML, platform engineering, ML compiler teams).
The project is specific enough to discuss in a technical interview, broad enough to connect to multiple domains, and concrete enough that you can actually demo it running.
Architecture Overview
The engine consists of four tightly coupled components. Here’s how they fit together:
┌─────────────────────────────────────────────────────┐
│ MODEL WEIGHTS │
│ (FP16/FP32 original, quantized offline) │
└──────────────────────┬──────────────────────────────┘
│ Packed into int4 nibbles
▼
┌─────────────────────────────────────────────────────┐
│ PACKED-NIBBLE WEIGHT STORE │
│ uint8 buffer: two int4 values per byte │
│ Paired with per-block scale & zero-point arrays │
└──────────────────────┬──────────────────────────────┘
│ During inference
▼
┌─────────────────────────────────────────────────────┐
│ FUSED DEQUANTIZING GEMM KERNEL │
│ 1. Read nibbles → unpack to int4 │
│ 2. Dequantize: (qval - zero_point) * scale │
│ 3. Accumulate into output via GEMM │
│ (single pass, no intermediate FP16 buffer) │
└──────────────────────┬──────────────────────────────┘
│
▼
┌─────────────────────────────────────────────────────┐
│ OUTPUT ACTIVATIONS │
│ (FP16 or BF16, ready for next layer or output) │
└─────────────────────────────────────────────────────┘
The key architectural decision is the fusion boundary: by merging dequantization and matrix multiplication into a single loop nest, we eliminate the intermediate tensor that would otherwise consume precious memory bandwidth. For a 7B parameter model quantized to INT4, the weight tensor is ~3.5GB. A naive dequantize-first approach temporarily doubles that to ~7GB during each GEMM. Fused execution avoids this entirely.
Component breakdown:
- Weight Packer: Takes original FP16 weights, computes per-group quantization parameters, and packs the resulting int4 values into a compact uint8 buffer.
- Quantization Config: Stores block size (typically 128 or 256), number of groups, and metadata for each weight tensor.
- Fused GEMM Kernel: The computational heart. Iterates over output rows and columns, dequantizes weights on-the-fly from nibble-packed storage, and accumulates partial sums.
- Inference Runtime: Orchestrates the forward pass across layers, manages memory, and provides a simple API for running a model.
Building It Step by Step
We’ll implement this in Python with NumPy for the reference logic, then show the critical kernel in C with SIMD intrinsics for the production-grade version. Python gets you running fast; C shows what actually matters for performance.
Step 1: Nibble Packing and Unpacking
Weights are stored two int4 values per byte. The high nibble holds the first value, the low nibble the second.
import numpy as np
from typing import Tuple
def pack_int4(values: np.ndarray) -> Tuple[np.ndarray, np.ndarray]:
"""
Pack a 1D array of int4 values (0..15) into a uint8 buffer.
Returns (packed_bytes, odd_mask) where odd_mask tracks
whether each original index was odd-positioned.
"""
assert values.dtype == np.int8 or values.dtype == np.uint8
assert values.max() <= 15 and values.min() >= 0
# Ensure even length
if len(values) % 2 != 0:
values = np.append(values, 0)
packed = np.zeros(len(values) // 2, dtype=np.uint8)
# High nibble = even index, low nibble = odd index
packed[:] = (values[0::2].astype(np.uint8) << 4) | values[1::2].astype(np.uint8)
return packed
def unpack_int4(packed: np.ndarray) -> np.ndarray:
"""Unpack a uint8 buffer back to int4 values as int8."""
high = ((packed >> 4) & 0x0F).astype(np.int8)
low = (packed & 0x0F).astype(np.int8)
# Interleave
result = np.empty(len(packed) * 2, dtype=np.int8)
result[0::2] = high
result[1::2] = low
return result
Step 2: Blockwise Quantization with Per-Group Scales and Zero-Points
We use the standard affine quantization formula: q = round(x / scale + zero_point), with q clamped to [0, 15]. The scale and zero-point are computed per group of block_size elements.
def quantize_blockwise(
weights: np.ndarray,
block_size: int = 128
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
"""
Quantize FP16 weights to INT4 blockwise.
Returns (packed_weights, scales, zero_points).
"""
flat = weights.reshape(-1)
num_groups = (len(flat) + block_size - 1) // block_size
scales = np.zeros(num_groups, dtype=np.float32)
zero_points = np.zeros(num_groups, dtype=np.int32)
for g in range(num_groups):
start = g * block_size
end = min(start + block_size, len(flat))
group = flat[start:end]
# Compute scale and zero-point per group
g_min = group.min()
g_max = group.max()
# Scale maps [g_min, g_max] to [0, 15]
scale = (g_max - g_min) / 15.0
if scale < 1e-8:
scale = 1e-8
# Zero-point: dequantized 0 maps to this value
zero_point = -g_min / scale
zero_point = np.clip(np.round(zero_point), 0, 15).astype(np.int32)
scales[g] = scale
zero_points[g] = zero_point
# Quantize
q = np.round(group / scale + zero_point).astype(np.int8)
q = np.clip(q, 0, 15)
flat[start:end] = q
packed = pack_int4(flat.astype(np.uint8))
return packed, scales, zero_points
Step 3: The Fused Dequantizing GEMM Kernel (C with SIMD)
This is the core of the project. The kernel reads nibbles directly from packed memory, dequantizes on the fly, and accumulates — never materializing the full FP16 weight matrix.
#include <stdint.h>
#include <immintrin.h> // AVX2 intrinsics
// Fused dequantizing GEMM: D = relu(A @ dequantize(Q))
// Q: packed int4 weights (uint8), scales and zero_points per group
// A: input activation matrix (M x K, row-major)
// D: output matrix (M x N, row-major)
// K is the reduced dimension; N is the output dimension
// block_size: elements per quantization group (e.g., 128)
void fused_dequant_gemm(
const uint8_t* __restrict__ Q_packed,
const float* __restrict__ scales,
const int32_t* __restrict__ zero_points,
const float* __restrict__ A,
float* __restrict__ D,
int M, int N, int K,
int block_size
) {
// K must be divisible by block_size for simplicity
int num_groups = K / block_size;
for (int m = 0; m < M; m++) {
for (int n = 0; n < N; n++) {
float acc = 0.0f;
const float* a_row = A + m * K;
for (int g = 0; g < num_groups; g++) {
float scale = scales[g];
int32_t zp = zero_points[g];
int g_start = g * block_size;
// Process block_size elements
// Unpack nibbles and dequantize inline
for (int k = 0; k < block_size; k += 2) {
uint8_t byte = Q_packed[g_start / 2 + k / 2];
int8_t q0 = (byte >> 4) & 0x0F;
int8_t q1 = byte & 0x0F;
// Dequantize: x = (q - zp) * scale
float dq0 = (float)(q0 - zp) * scale;
float dq1 = (float)(q1 - zp) * scale;
acc += a_row[g_start + k] * dq0;
acc += a_row[g_start + k + 1] * dq1;
}
}
// ReLU activation
D[m * N + n] = acc > 0.0f ? acc : 0.0f;
}
}
}
For the AVX2-optimized inner loop, you’d vectorize the unpack-dequantize-accumulate pattern using 256-bit registers. The key insight is that you can load 16 packed bytes (32 int4 values), unpack them using vpsrlvd and vpand, then multiply-add against 16 activation values in a single pass.
// AVX2-optimized inner block (conceptual sketch)
// Processes 16 output channels simultaneously
__m256 compute_block_avx2(
const uint8_t* q_ptr,
const float* a_ptr,
float scale,
int32_t zp
) {
__m256i q_bytes = _mm256_loadu_si256((__m256i*)q_ptr);
// Unpack high nibbles: shift right by 4, mask with 0x0F
__m256i q_high = _mm256_and_si256(_mm256_srli_epi16(q_bytes, 4), _mm256_set1_epi16(0x0F));
__m256i q_low = _mm256_and_si256(q_bytes, _mm256_set1_epi16(0x0F));
// Convert to float and dequantize
__m256 q_high_f = _mm256_cvtepi32_ps(_mm256_cvtepu16_epi32(_mm256_packus_epi16(q_high, q_low)));
// ... (full unpacking to float)
__m256 scale_v = _mm256_set1_ps(scale);
__m256 zp_v = _mm256_set1_ps((float)zp);
__m256 dq = _mm256_sub_ps(q_high_f, zp_v);
dq = _mm256_mul_ps(dq, scale_v);
// Load activations and multiply-accumulate
__m256 a_v = _mm256_loadu_ps(a_ptr);
return _mm256_fmadd_ps(a_v, dq, _mm256_setzero_ps());
}
Step 4: Python Reference Inference Runtime
Wire everything together into a runnable inference loop:
class BlockwiseInt4Linear:
"""A single linear layer with INT4 blockwise weights."""
def __init__(self, in_features: int, out_features: int, block_size: int = 128):
self.in_features = in_features
self.out_features = out_features
self.block_size = block_size
self.num_groups = in_features // block_size
# Allocate weight storage
self.weight_packed = np.zeros(self.num_groups * block_size // 2, dtype=np.uint8)
self.scales = np.zeros(self.num_groups, dtype=np.float32)
self.zero_points = np.zeros(self.num_groups, dtype=np.int32)
self.bias = np.zeros(out_features, dtype=np.float32)
def load_quantized(self, fp16_weights: np.ndarray):
"""Quantize and store FP16 weights."""
packed, scales, zero_points = quantize_blockwise(fp16_weights, self.block_size)
self.weight_packed = packed
self.scales = scales
self.zero_points = zero_points
def forward(self, x: np.ndarray) -> np.ndarray:
"""Forward pass: y = x @ W^T + b, with INT4 weights."""
# Use the fused GEMM logic (Python reference)
K = self.in_features
M = x.shape[0]
N = self.out_features
# Reconstruct the dequantized weight matrix on-the-fly
W_dq = np.zeros((N, K), dtype=np.float32)
flat_q = unpack_int4(self.weight_packed)
for g in range(self.num_groups):
start = g * self.block_size
end = start + self.block_size
scale = self.scales[g]
zp = self.zero_points[g]
W_dq[:, start:end] = (flat_q[start:end].astype(np.float32) - zp) * scale
# Add bias and apply activation
output = x @ W_dq.T + self.bias
return np.maximum(output, 0.0) # ReLU
# --- End-to-end test ---
if __name__ == "__main__":
# Create a simple 2-layer model
layer1 = BlockwiseInt4Linear(in_features=256, out_features=128)
layer2 = BlockwiseInt4Linear(in_features=128, out_features=10)
# Random FP16 weights
np.random.seed(42)
w1 = np.random.randn(128, 256).astype(np.float16)
w2 = np.random.randn(10, 128).astype(np.float16)
layer1.load_quantized(w1)
layer2.load_quantized(w2)
# Run inference
x = np.random.randn(4, 256).astype(np.float32)
h = layer1.forward(x)
y = layer2.forward(h)
print(f"Output shape: {y.shape}")
print(f"Output values:\n{y}")
Step 5: Build the C Extension
For the production path, compile the C kernel as a Python extension using cffi or pybind11:
# setup.py for pybind11 binding
# pip install pybind11 numpy
// bindings.cpp
#include <pybind11/pybind11.h>
#include <pybind11/numpy.h>
#include <cstdint>
namespace py = pybind11;
void fused_dequant_gemm(
const uint8_t* Q_packed,
const float* scales,
const int32_t* zero_points,
const float* A,
float* D,
int M, int N, int K,
int block_size
);
PYBIND11_MODULE(int4_engine, m) {
m.def("fused_gemm", &fused_dequant_gemm,
"Fused dequantizing GEMM with INT4 blockwise weights",
py::arg("Q_packed"), py::arg("scales"), py::arg("zero_points"),
py::arg("A"), py::arg("D"),
py::arg("M"), py::arg("N"), py::arg("K"),
py::arg("block_size"));
}
# Compile
c++ -O3 -march=native -O3 -shared -std=c++17 -fPIC \
$(python3 -m pybind11 --includes) bindings.cpp -o int4_engine$(python3-config --extension-suffix)
Running and Testing It
Local Setup
# Clone and set up
git init int4-engine
cd int4-engine
python3 -m venv venv && source venv/bin/activate
pip install numpy pybind11
Unit Tests
Verify correctness against a full-precision reference:
# test_engine.py
import numpy as np
from int4_engine import fused_gemm
def test_correctness():
"""Compare fused INT4 GEMM output against FP32 reference."""
np.random.seed(123)
M, N, K = 4, 8, 256
block_size = 128
num_groups = K // block_size
# Create reference weights and activations
W_fp32 = np.random.randn(N, K).astype(np.float32) * 0.5
A = np.random.randn(M, K).astype(np.float32)
b = np.random.randn(N).astype(np.float32)
# Quantize weights
packed, scales, zero_points = quantize_blockwise(W_fp32.flatten(), block_size)
W_int4 = W_fp32.copy()
flat_q = unpack_int4(packed)
for g in range(num_groups):
start = g * block_size
end = start + block_size
W_int4[:, start:end] = (flat_q[start:end].astype(np.float32) - zero_points[g]) * scales[g]
# Reference output
ref_output = A @ W_int4.T + b
ref_output = np.maximum(ref_output, 0.0)
# Fused kernel output
D = np.zeros((M, N), dtype=np.float32)
fused_gemm(
packed, scales, zero_points,
A, D,
M, N, K, block_size
)
# Check accuracy
max_err = np.max(np.abs(D - ref_output))
mean_err = np.mean(np.abs(D - ref_output))
print(f"Max absolute error: {max_err:.6f}")
print(f"Mean absolute error: {mean_err:.6f}")
assert max_err < 1e-3, f"Accuracy too low: {max_err}"
print("✅ Correctness test passed!")
def test_throughput():
"""Benchmark the fused kernel against a naive dequantize-then-GEMM."""
import time
M, N, K = 64, 512, 4096
block_size = 128
num_groups = K // block_size
A = np.random.randn(M, K).astype(np.float32)
packed, scales, zero_points = quantize_blockwise(
np.random.randn(N, K).astype(np.float32) * 0.5, block_size
)
# Fused kernel
D_fused = np.zeros((M, N), dtype=np.float32)
start = time.perf_counter()
fused_gemm(packed, scales, zero_points, A, D_fused, M, N, K, block_size)
fused_time = time.perf_counter() - start
# Naive: dequantize then GEMM
W_dq = np.zeros((N, K), dtype=np.float32)
flat_q = unpack_int4(packed)
for g in range(num_groups):
s = g * block_size
W_dq[:, s:s+block_size] = (flat_q[s:s+block_size].astype(np.float32) - zero_points[g]) * scales[g]
start = time.perf_counter()
D_naive = A @ W_dq.T
naive_time = time.perf_counter() - start
print(f"Fused GEMM: {fused_time*1000:.2f}ms")
print(f"Naive dequantize+GEMM: {naive_time*1000:.2f}ms")
print(f"Speedup: {naive_time/fused_time:.2f}x")
print(f"Memory savings: avoided allocating {N*K*4/1e6:.1f}MB FP32 weight buffer")
if __name__ == "__main__":
test_correctness()
test_throughput()
python test_engine.py
Expected output shows the fused kernel matching the reference within floating-point tolerance while avoiding the intermediate weight allocation entirely. The throughput test should show meaningful speedup on larger matrices, especially as the activation batch size grows.
Extending It: Your Roadmap to Senior-Level
The baseline engine above is a strong portfolio piece. But to signal senior-level capability, layer on these production-grade upgrades:
Model Persistence and Checkpointing — Implement a binary serialization format (think: a simplified version of GGML’s file format) that stores quantized weights, quantization metadata, and architecture config in a single
.int4file. This matters because production systems need reproducible deploys and versioned model artifacts. Write asave()/load()API that handles endianness, checksums, and backward compatibility.Batched Inference and Horizontal Scaling — Add a gRPC or HTTP service layer (using something like Triton Inference Server as inspiration) that accepts batched requests, schedules them across multiple worker threads with thread-local weight buffers, and returns results with latency headers. This signals you understand how inference engines actually serve traffic in production, not just run locally.
Observability and Profiling — Integrate structured logging and metrics collection (Prometheus client, OpenTelemetry spans) that tracks per-request latency percentiles, memory bandwidth utilization, and quantization error per layer. A senior engineer doesn’t just build the engine — they build the instrumentation that tells you why it’s slow or wrong.
Fault Tolerance and Graceful Degradation — Add fallback paths: if INT8 or FP16 kernels are available on the host hardware (detected via CPUID), automatically use them for layers where INT4 quantization error exceeds a configurable threshold. This demonstrates systems thinking about reliability — the engine should never silently produce garbage output.
Kernel Auto-Tuning and Benchmarking Harness — Build a benchmarking framework that auto-tunes kernel parameters (block size, tile dimensions, vector width) against the specific CPU microarchitecture at deploy time. Store results in a cache so subsequent runs skip the search. This is the kind of work that appears in ML compiler teams at companies like Meta and Google.
CUDA/HIP Backend Port — Extend the fused kernel to GPU using CUDA or HIP, leveraging shared memory for activation tiling and cooperative fetching for weight nibbles. This is the single biggest signal upgrade — it proves you can navigate the GPU programming model, understand memory coalescing, and write kernels that actually run fast on NVIDIA hardware.
Each of these upgrades maps directly to a real production concern: deployment (persistence), scalability (batching), operational confidence (observability), reliability (fault tolerance), performance engineering (auto-tuning), and hardware acceleration (GPU).
Key Takeaways
- Packed-nibble storage cuts model weight size by 2x compared to INT8 and 4x compared to FP16, and the dequantize-and-multiply pattern is the fundamental building block of every modern LLM inference engine.
- Fusing dequantization into GEMM eliminates the intermediate FP16/FP32 weight buffer, which is often the binding constraint on memory bandwidth — the real bottleneck in inference, not compute.
- Per-group scales and zero-points are essential because a single global quantization parameter cannot capture the dynamic range of all weight channels; blockwise parameters keep quantization error bounded.
- This project sits at the exact intersection of ML and systems engineering — the same skills are needed for building inference engines, ML compilers, database query optimizers, and embedded accelerators.
- The code you write here is directly portable to CUDA, Metal, or WebGPU backends, and the quantization logic is identical to what frameworks like llama.cpp and GGUF use in production.
- Each extension (persistence, scaling, observability, GPU) adds a dimension to your portfolio that signals readiness for real production systems, not just academic prototypes.
Further Reading
- LLM.int8(): 8-bit Matrix Multiplication for Transformers at Scale — The foundational paper on mixed-precision decomposition that makes blockwise quantization work for large attention layers. Essential reading for understanding why per-group quantization is necessary.
- GGML and GGUF: The Quantization Format Behind llama.cpp — The canonical production implementation of blockwise INT4/INT8 quantization. Study their
ggml-quants.cfor the exact nibble-packing and dequantization kernels this project is based on. - CUDA Programming Guide: Memory Hierarchy and Coalesced Access — When you port the fused kernel to GPU, this is the reference for understanding shared memory, global memory coalescing, and the memory bandwidth constraints that dominate inference latency.
- Block-Sparse and Blockwise Quantization for Neural Networks (ICLR 2023) — A rigorous treatment of why blockwise (as opposed to per-tensor or per-channel) quantization achieves the best accuracy-efficiency tradeoff, with theoretical bounds on quantization error.
- [The Matrix Multiplication Algorithm Handbook](https://www.cs.utexas.edu/~flame/pubs/goto paper.pdf) — The Goto paper on an optimized GEMM implementation. Understanding its tiling, blocking, and vectorization strategies will directly inform how you structure the fused dequantize-GEMM loop nest.
- OpenTelemetry Specification — When you add observability to your engine, this is the standard for instrumenting metrics, traces, and logs in a way that integrates with modern observability stacks.
This project gives you something rare in a portfolio: a single codebase that demonstrates depth across numerical methods, memory architecture, kernel optimization, and production systems thinking. Start with the Python reference, verify correctness, then push the C kernel and add one extension at a time. Each layer you add is a concrete conversation starter in an interview — and a testament to engineering skill that hiring managers can evaluate immediately.