LAB / GPU PROGRAMMING

DAY 15 / 30 · WRITING CUDA · ~12 MIN READ

Tiling with shared memory

Speeding up matmul by staging tiles in shared memory instead of hitting global memory per element.

01 · The idea in one sentence

Load a tile once, use it many times

Day 12's verdict on naive matmul: 549 GB of DRAM traffic for 137 GFLOPs of math, intensity 0.25 FLOP/byte, ceiling 0.4% of peak — because every input element is re-fetched N times. Day 08 built the tool that fixes this: shared memory, a per-block scratchpad. The fix: have each block compute a tile of the output, and stage the matching tiles of A and B in shared memory so every loaded byte feeds all 256 of the block's threads.

The block loops over the K dimension in phases. In each phase it cooperatively loads one 16×16 tile of A and one of B — coalesced, since consecutive threads take consecutive addresses — syncs, computes 16 multiply-accumulates per output element from SRAM, and moves to the next phase. Nothing about the math changes; only the data flow does.

02 · The tiled kernel

Twenty lines that change the roofline

#define T 16

__global__ void matmul_tiled(const float* A, const float* B,
                                   float* C, int N){
    // note the padding from day 08 — column reads stay conflict-free
    __shared__ float As[T][T + 1];
    __shared__ float Bs[T][T + 1];

    int tx = threadIdx.x, ty = threadIdx.y;
    int row = blockIdx.y * T + ty;  // this thread's output element
    int col = blockIdx.x * T + tx;
    float acc = 0.0f;

    // loop over K in phases of T
    for (int k0 = 0; k0 < N; k0 += T){
        // each thread loads ONE element of each tile → coalesced
        As[ty][tx] = A[row * N + (k0 + tx)];
        Bs[ty][tx] = B[(k0 + ty) * N + col];
        __syncthreads();  // tiles fully staged

        // 16 MACs per thread, all served from SRAM
        for (int k = 0; k < T; ++k)
            acc += As[ty][k] * Bs[k][tx];

        __syncthreads();  // done reading before next phase overwrites
    }
    C[row * N + col] = acc;
}

Three details carry the correctness. The first __syncthreads() is day 08's fence — no thread reads a tile until all 256 have written theirs. The second sync is easier to forget: it stops a fast thread from overwriting the tiles while a slow thread is still reading the previous phase. And the +1 padding is day 08's bank-conflict fix — the inner loop's As[ty][k] reads a row (fine) but Bs[k][tx] reads down a column, which without the pad would serialize 16-way. (This version assumes N is a multiple of 16; production kernels add day 14's guards and a tiled epilogue.)

03 · The arithmetic of what we saved

16× less traffic, same math

Count traffic per block for one full sweep of K. Naive: each of 256 output elements reads 2N floats → 512N floats per block. Tiled: each phase loads 2 × 16 × 16 floats, and there are N/16 phases → 2 × 256 × N/16 = 32N floats. Same outputs, 16× less DRAM traffic. Across the whole problem, 549 GB drops to ~34 GB, and day 11's chart moves: intensity 0.25 → 4 FLOP/byte, ceiling 126 GFLOPS → ~2,000 GFLOPS.

NAIVE · TRAFFIC

549 GB → ceiling ~126 GFLOPS (0.4%)

TILED 16×16 · TRAFFIC

~34 GB → ceiling ~2,016 GFLOPS (6.9%)

TILED 32×32 · TRAFFIC

~17 GB → ceiling ~8,064 GFLOPS (28%)

THE PATTERN

reuse per byte = tile size: T× the naive ceiling

On real hardware the measured speedup for this exact rewrite is typically 8–16× — the same order as the traffic math, minus overheads the model doesn't count. Not the full 16× in wall-clock, because the kernel is still memory-bound and still pays for the phases' synchronization — but the direction and magnitude are day 11's roofline made visible.

04 · Where this plateaus

Why the pros keep climbing and you shouldn't (yet)

The table shows the ladder: bigger tiles, more reuse — 32×32 tiles quadruple the ceiling again. But each rung costs resources day 06 and day 09 priced: a 32×32 tile needs 8 KB of shared memory per stage, and each thread accumulates 32 registers' worth of outputs before epilogue. Push tile size and register pressure until occupancy craters, and the gains evaporate. Production GEMM libraries layer on double buffering (load phase k+1 while computing k), register tiling (each thread computes a small sub-tile), and tensor-core instructions — reaching 80–90% of peak where our 16×16 tile sits at a few percent.

Don't take that as defeat. You now understand the mechanism every one of those tricks exploits: keep bytes on-chip longer. The method — count traffic, tile for reuse, watch the registers — is exactly what you'll use to read those libraries on day 27, and exactly what day 24's FlashAttention does for attention.

CHEATSHEET · DAY 15
  • 01Tile the output; stage the inputs in shared memory. Each block computes a T×T output tile by sweeping K in phases, loading A- and B-tiles once per phase. Reuse per byte becomes ~T.
  • 02Two __syncthreads() per phase, both required. One after staging (don't read unfinished tiles), one after computing (don't overwrite tiles still being read).
  • 03Pad shared tiles: [T][T+1]. Bs[k][tx] walks columns; without padding the 16-wide stride lands on one bank and serializes (day 08).
  • 04Tile size is a resource negotiation. Bigger tiles mean more reuse but more shared memory and registers per block — the day 06/09 budgets decide where it stops paying.

Next: day 16 — parallel reductions. Summing a million numbers looks trivial and is quietly one of the most instructive kernels you'll ever write.

get in touch