GPU Optimization
Where AI performance actually comes from — occupancy, memory access patterns, fusion, and precision.
Before touching code, classify the kernel: is it memory-bound (limited by HBM bandwidth) or compute-bound (limited by ALU/tensor-core throughput)? The roofline model further down answers this precisely, but a rough gut-check is usually enough to start — if a kernel does little arithmetic per byte loaded (elementwise ops, normalization, most attention kernels at long sequence length), assume memory-bound; if it does a lot of reuse per byte (large GEMMs, convolutions), assume compute-bound. Optimizing the wrong side of that line is the single most common way to spend an afternoon and gain nothing: fusing kernels doesn't help a compute-bound GEMM, and switching to lower precision doesn't help a kernel that was never bandwidth-limited to begin with. The sections below are ordered from foundational (works on almost any kernel) to advanced (situational, but where the largest remaining wins usually live once the basics are done) — profile with the tools in the last section before and after every change, since intuition about GPUs is wrong often enough to be worth distrusting.
Occupancy is the ratio of active warps/wavefronts per SM/CU to the hardware maximum. High occupancy is how a GPU hides memory latency — while one warp waits on a load, another warp gets scheduled (see the warp scheduler). Occupancy is capped by whichever resource runs out first per block: registers, shared memory, or thread-block slots. A kernel that uses too many registers per thread silently reduces the number of warps that can be resident, which lowers occupancy and can hurt latency hiding — even though the extra registers were "helping" that one thread.
All threads in a warp/wavefront execute the same instruction each cycle. When an if/else takes different branches within the same warp, the hardware serializes: it runs the true-branch threads while masking off the false-branch threads, then vice versa — so divergent branches cost the sum of both paths, not the max. This is why GPU code favors branchless formulations (predication, min/max, masks) over conditionals inside hot loops whenever the branch outcome varies across neighboring threads.
When consecutive threads in a warp access consecutive addresses in global memory, the hardware merges (coalesces) them into one wide memory transaction. When accesses are scattered (e.g. a strided or transposed access pattern), the same data volume costs many more transactions. This single access-pattern detail is one of the largest, most common sources of "why is my kernel memory-bound and slow" — often a bigger effect than any compute optimization.
Shared memory (Nvidia) / LDS (AMD) is on-chip, programmer-managed, and roughly two orders of magnitude faster than global memory — but it's split into banks, and if multiple threads in the same warp hit the same bank on different addresses simultaneously, those accesses serialize (a bank conflict). Classic mitigation: pad a shared-memory array's row stride by one element to break the power-of-two access pattern that causes conflicts in matrix-transpose-style kernels.
Every unfused op (e.g. matmul, then bias-add, then activation) is a separate kernel launch that round-trips its output through global memory before the next kernel reads it back. Fusion combines multiple ops into one kernel so intermediate results stay in registers/shared memory, cutting both launch overhead and memory traffic. This is exactly what torch.compile's Inductor backend and hand-written epilogues in cuBLASLt/hipBLASLt are doing — folding bias/activation/scaling into the GEMM kernel itself.
| Format | Exponent/mantissa bits | Trade-off |
|---|---|---|
| FP32 | 8 / 23 | full range and precision, slowest, most memory |
| FP16 | 5 / 10 | narrow dynamic range — prone to overflow/underflow, needs loss scaling |
| BF16 | 8 / 7 | same range as FP32 (no loss scaling needed), less mantissa precision — the default for most modern training |
| FP8 | varies (e4m3 / e5m2) | half the memory/bandwidth of FP16 again, used increasingly for inference and some training with per-tensor scaling |
Lower precision helps in two places at once: it doubles (or more) the FLOPs a tensor/matrix core can issue per cycle, and it halves the memory traffic to move the same tensor — which matters most since most AI kernels are memory-bound (see architecture).
Given a kernel's arithmetic intensity (FLOPs/byte), the roofline model tells you which optimization is worth attempting. Below the ridge point: reduce memory traffic (fusion, better access patterns, lower precision) — adding more compute won't move the needle. Above the ridge point: better tiling/instruction scheduling and higher-throughput math (tensor cores, lower precision) is what helps. Profiling without first checking which regime a kernel is in is the most common way to optimize the wrong thing.
The classic advanced technique, best seen in matmul: instead of every thread re-reading its row/column of the input matrices from global memory for every output element, threads in a block cooperatively load small tiles of each matrix into shared memory once, then every thread in the block reuses those tiles for many multiply-accumulates before the block moves on to the next tile. This trades global-memory bandwidth for shared-memory bandwidth (roughly 100× faster) and is the single biggest reason a naive matmul kernel and a well-tuned one can differ by an order of magnitude, even though both do the exact same number of FLOPs. cuBLAS/rocBLAS, CUTLASS, and Triton all pick tile sizes as their central tuning knob for exactly this reason — too small and you don't amortize the global-memory load, too large and you run out of shared memory or registers and occupancy drops.
Occupancy hides latency within a kernel by overlapping warps; streams hide latency across operations by overlapping a kernel with a memory copy, or one kernel with another, on independent hardware queues. The standard pattern is double buffering: while the GPU computes on chunk N, an async copy (see CUDA & HIP) is already bringing chunk N+1 in on a separate stream, so compute never stalls waiting on PCIe or host-to-device transfer. This matters most for pipelines that can't fit the whole problem in VRAM at once (large dataset streaming, KV-cache paging) or that mix meaningfully sized host/device transfers with compute, rather than for a single self-contained kernel launch where there's nothing else to overlap with.
Past a certain point, hand-picking tile sizes, unroll factors, and pipeline stages by intuition stops working — the optimal choice depends on the exact problem shape, GPU generation, and precision, and shifts between all of them. Frameworks like Triton and CUTLASS lean on autotuning instead: define the kernel as a template over these parameters, benchmark a search space of concrete configurations on the real hardware, and cache whichever one wins for a given input shape. This is why the first call to a Triton kernel with a new shape is often slow (it's compiling and benchmarking candidates) while subsequent calls hit a cached, tuned version — the same idea behind PyTorch's torch.backends.cudnn.benchmark = True for convolution algorithm selection.