TL;DR: A naive GEMM has one thread per output element looping over K, loading one element of A and one of B from global memory per multiply-add: intensity of 0.25 FLOP/B against a ridge near 300, so it runs at a few hundred GFLOPS. Shared-memory tiling loads a 128 × 32 tile of A and 32 × 128 of B once per K-step and lets 128 × 128 outputs reuse them, raising intensity to about 64 FLOP/B and speed by 20 to 30×. Register blocking (each thread computes an 8 × 8 micro-tile from registers) cuts shared-memory traffic by another 8× and reaches about 15 to 20 TFLOPS in fp32, near the fp32 CUDA-core peak. Vectorized 16-byte loads and double buffering (prefetch the next tile while computing this one, with
cp.asyncor TMA) push toward the CUDA-core ceiling. Tensor cores change the roof:mmaorwgmmainstructions with the tiles in the layouts they require reach 500 to 700 TFLOPS in bf16, 50 to 70% of the 990 peak, with the remaining gap in pipelining depth, warp specialization and epilogue overlap that CUTLASS and cuBLAS encode.
How to approach it
State the ridge and the naive kernel's intensity to place it. Then each step: what it adds, the new intensity or the new bottleneck, and a realistic speedup. Show the register-blocked tile's inner loop. Then tensor cores as the change of roof. Then the last 30% and where to stop.
A strong answer
A typical situation: an interview take-home asks for a bf16 GEMM at 4,096³ and a report of TFLOPS at each optimization; the candidate who shows the intensity arithmetic beside each measurement is the one the reviewers remember.
The reference numbers on an H100 SXM:
peak: bf16 tensor cores ~990 TFLOPS dense; fp32 CUDA cores ~67 TFLOPS; HBM ~3.35 TB/s
ridge: 990e12 / 3.35e12 ≈ 295 FLOP/B (tensor cores); 67e12 / 3.35e12 ≈ 20 FLOP/B (fp32 cores)
the problem: M = N = K = 4,096; FLOPs = 2 × 4,096³ ≈ 1.37e11; A and B in bf16 = 32 MB each, C = 32 MB
minimum HBM traffic ≈ 96 MB → the whole GEMM's intensity from HBM is 1,430 FLOP/B, far above the ridge; the GEMM is
compute-bound in principle; every step below is about making the on-chip traffic match that
Step by step:
0. naive: thread (i, j) loops k = 0..K−1: acc += A[i][k] × B[k][j]
per MAC: 2 loads from global (or L1/L2) of 2 bytes each, 2 FLOPs → 0.5 FLOP/B; the B access is strided by N
(uncoalesced) → effective intensity lower still; result: ~0.5 to 2 TFLOPS in fp32 (a few percent of even the
CUDA-core peak); it is a bandwidth- and latency-bound kernel
1. shared-memory tiling (block tile BM × BN = 128 × 128, K-step BK = 32)
each K-step: load A tile 128 × 32 and B tile 32 × 128 (16 KB in bf16) cooperatively and coalesced into shared
memory; compute 128 × 128 × 32 MACs from shared
global intensity: 2 × 128 × 128 × 32 FLOPs / 16 KB ≈ 64 FLOP/B → above the fp32 ridge; the kernel is now bounded
by shared-memory bandwidth and instruction issue, not HBM; measured ~8 to 12 TFLOPS fp32 (20 to 30× over naive)
2. register blocking (thread tile TM × TN = 8 × 8, 256 threads per block)
each thread keeps an 8 × 8 accumulator in 64 registers; per k it loads 8 values of A and 8 of B from shared into
registers and does 64 MACs → shared traffic per MAC drops 8× (from 2 loads per MAC to 16 loads per 64 MACs);
shared-memory bandwidth is no longer the limit; measured ~15 to 20 TFLOPS fp32, near the CUDA-core ceiling
this is the step that requires the bank-conflict-free layout: the A tile is stored transposed (or padded) so the
8 values a thread needs are in distinct banks
3. vectorized loads and stores (16-byte float4 / 8 × bf16 per load), and an output epilogue that writes C in
coalesced 16-byte chunks; fewer instructions, closer to the ceiling; +10 to 20%
4. double buffering (two shared-memory tile buffers; issue cp.async for tile k+1 before computing tile k; wait and
swap): hides global-memory latency behind compute; on Hopper, TMA (a single instruction that moves a whole tile
with the swizzle applied) does this with less register and instruction cost; +10 to 30% depending on how exposed
the latency was; with steps 1 to 4 a careful fp32 kernel reaches ~60 to 70% of the fp32 CUDA-core peak
5. tensor cores: replace the inner FMAs with mma.sync (Ampere-style, warp-level 16 × 8 × 16 tiles) or wgmma (Hopper,
warpgroup-level 64 × N × 16 with operands in shared memory in a swizzled layout); the accumulator stays in
registers; the roof rises from 67 to 990 TFLOPS; the kernel's job becomes feeding the tensor cores: the
shared-memory tile layout must match what ldmatrix / wgmma expect, and the K-loop must keep several tiles in
flight; measured: a first working wgmma kernel ~300 to 400 TFLOPS; with TMA, a 3 to 4 stage pipeline and a
fused epilogue ~500 to 700 TFLOPS (50 to 70% of peak)
sanity: the intensity of a 128 × 128 × 64 block tile in bf16 is 2 × 128 × 128 × 64 / (2 × 128 × 64 × 2 B) ≈ 64
FLOP/B from shared memory per tile; the tensor cores at 990 TFLOPS need ~15 TB/s from shared, and an SM's
shared memory delivers on the order of that only with wgmma reading operands directly, which is why the
Hopper path is designed the way it is.
Tiled Matrix Multiplication has the tile hierarchy; Shared Memory and Bank Conflicts is why step 2's layout matters; Tensor Cores and Matrix Units is what step 5 targets; Roofline Model is the frame for every number above.
The register-blocked inner loop, in outline:
// per K-step (BK = 32), after the block's A tile (BM×BK, stored transposed: As[BK][BM]) and B tile (Bs[BK][BN]) are in shared
float acc[8][8] = {0};
for (int k = 0; k < BK; ++k) {
float a[8], b[8];
for (int i = 0; i < 8; ++i) a[i] = As[k][ty * 8 + i]; // 8 consecutive floats: one float4 pair, conflict-free
for (int j = 0; j < 8; ++j) b[j] = Bs[k][tx * 8 + j];
for (int i = 0; i < 8; ++i)
for (int j = 0; j < 8; ++j)
acc[i][j] += a[i] * b[j]; // 64 FMAs per 16 shared loads
}
The last 30%, and where to stop: past 70% the gains come from pipelining depth (4 or more stages so the tensor cores never wait), warp specialization (producer warps issue TMA loads while consumer warpgroups issue wgmma), persistent kernels that loop over output tiles to avoid launch and tail effects, split-K for small M and N, and epilogue overlap; each is a few percent and each is a page of code. That is what CUTLASS templates and cuBLAS heuristics contain (CUTLASS and Tensor Core Kernels), and the interview answer for production is "use them, and write the kernel yourself when you need an epilogue or a layout they do not offer".
The reversal condition: for small or skinny GEMMs (M = 8 at decode time), the kernel is bandwidth-bound on the weight read regardless of tiling, and the optimization target is reading the weights once at full bandwidth, which is a different kernel (a batched GEMV) with none of the tiling above. Nsight Compute's tensor-pipe and memory rows are how you tell which of the six steps is still leaving performance behind.
What interviewers probe next
- "Why 128 × 128 and not 256 × 256?" Larger tiles raise intensity but need more shared memory and registers per block, cutting occupancy and the ability to hide latency; the sweet spot is measured per architecture; Hopper kernels use 128 × 256 or 256 × 128 with wgmma.
- "What does the intensity of the whole GEMM tell you?" That HBM is not the limit for large square GEMMs; the on-chip hierarchy is, which is why every step is about shared memory and registers.
- "How do you verify correctness?" Against cuBLAS with a tolerance scaled by K and dtype; test non-multiples of the tile size; check the accumulator is fp32.
- "What changes for bf16 inputs with fp32 accumulation on CUDA cores?" Convert on load; the FMAs are fp32; the traffic halves; the ceiling is still 67 TFLOPS until tensor cores are used.
Common mistakes
- Tiling without register blocking, and a kernel bound by shared-memory bandwidth at 10 TFLOPS.
- Claiming a tensor-core kernel is "at peak" at 300 TFLOPS on an H100.
- A tile layout with bank conflicts, and the register-blocking step giving nothing.
- Measuring a single warm run and reporting it as throughput.
Key takeaways
- Naive: 0.5 FLOP/B and ~1% of peak. Tiling: 64 FLOP/B from global. Register blocking: 8× less shared traffic. Vectorization and double buffering: to the CUDA-core ceiling.
- Tensor cores raise the roof 15×; the kernel becomes a feeding problem: layouts, TMA, pipeline stages.
- 50 to 70% of 990 TFLOPS is a strong hand-written bf16 kernel; the rest is warp specialization and persistence, which live in CUTLASS.
