Why Fast GPUs Still Wait for Memory

One operation, two costs

The roofline is a bound, and a bound is only as useful as the assumptions underneath it are visible. Four kinds of callout run through this article, and they are worth attending to rather than skipping — most of the ways this model gets misapplied in practice are written in them.

Nuance

A subtlety that changes how a result should be read — a bound that is tighter or looser than it looks, or a verdict that only holds at one boundary.

Trap

A conclusion the arithmetic invites but does not support. Each one is a mistake that is common enough to be worth naming.

Approximation

A number that is estimated, vendor-quoted or policy-dependent rather than derived. Where a decision turns on one of these, measure it yourself.

Aside

A worked example, a sanity check, or a practical suggestion. Skippable without losing the argument.

Take the simplest operation that does anything at all: add two vectors of n elements.

Elementwise addC = A + B

Count the arithmetic. One addition per element, so n floating-point operations. Now count the data. Two operands have to arrive and one result has to leave, so at s bytes per element the operation moves 3sn bytes. In bf16 that is 6n bytes to perform n operations — six bytes for every single FLOP.

Now do the same for a square matrix multiply of the same element type.

Matrix multiplyC = A·B, all matrices n × n

The arithmetic is 2n³: for each of the n² output elements, a dot product of length n, which is n multiplies and n adds. The data is three matrices of n² elements, so 3sn² bytes. In bf16 that is 6n² bytes to perform 2n³ operations — three bytes for every n FLOPs.

Both are "one operation" in the sense that a programmer writes one line. They are not remotely the same kind of thing. Put the ratios side by side:

FLOPs per byte elementwise add = n / 3sn = 1/(3s) = 0.17 (bf16, and independent of n)
matrix multiply = 2n³ / 3sn² = 2n/(3s) = n/3 (bf16, and it grows with n)

At n = 4096 those two numbers differ by a factor of eight thousand. No amount of counting operations will ever reveal that, because the operation counts are not what separates them. What separates them is how much arithmetic each does per byte it touches — and that quantity is the subject of everything that follows.

This matters because the machine has two independent capacities, and it publishes both. An H100 SXM5 will retire bf16 tensor-core arithmetic at roughly 989 TFLOP/s and stream data from HBM3 at roughly 3.35 TB/s. Two rates, two units, two separate parts of the chip that run at the same time. Every operation draws on both, in a ratio it does not get to choose, and finishes when the more heavily loaded of the two is done.

So the question this article answers is not "how many FLOPs does this cost" but: given an operation and a GPU, which of the two capacities runs out first? Answer that and you know what to fix. Get it wrong and you will spend a quarter optimising something whose time was never in your hands.

Approximation · both published rates are peaks, and one is a convention
989 TFLOP/s assumes every issue slot goes to a tensor-core matrix instruction at the boost clock; 3.35 TB/s assumes perfectly coalesced streaming. Neither holds by default, and the next section is about how far short of them a real kernel can fall. The FLOP count carries a convention too: a fused multiply-add is two operations here and throughout, which is the near-universal choice, but vendor headline figures sometimes assume 2:4 structured sparsity and are then exactly twice the dense number. Every peak quoted in this article is the dense one — if a comparison ever looks off by precisely two, check this first.

Figure 01

Stage 01 / 08

Figure 01 — Two pipes, one kernel. The memory pipe delivers bytes at β; the arithmetic array consumes operands at π. A kernel supplying I FLOPs per delivered byte can saturate at most one of them. Below the machine balance I_* = π/β the tensor cores bubble; above it the HBM link goes quiet.

Compute time, and why peak lies

Give the arithmetic a symbol. Let W be the number of floating-point operations an operation requires, and π the rate at which the machine can retire them. Then the time the arithmetic alone would take is a division:

Compute timeT_compute = W / π

That is the whole derivation, and it would be unremarkable except for one thing: π is not a number the machine owes you. It is the largest rate anyone has ever measured under conditions you are almost certainly not reproducing. The subtitle of this article is a claim about π, and it is worth making the claim precise before building anything on top of it.

π is the top of a stack

The published peak assumes every issue slot goes to a tensor-core matrix-multiply-accumulate instruction at the boost clock. Change the instruction mix and you are no longer measured against that line at all, but against a lower one:

Compute ceilings, H100 SXM5 fp8 tensor core     1979 TFLOP/s
bf16 tensor core    989 TFLOP/s   (2.0× below fp8)
tf32 tensor core     495 TFLOP/s
fp32 CUDA cores     67 TFLOP/s   (14.8× below bf16)
fp32 without FMA    33 TFLOP/s
transcendental (SFU) ≈17 TFLOP/s   (~¼ the fp32 rate)

A kernel that never issues an MMA cannot exceed the fp32 line however much arithmetic it has, and a kernel dominated by expf cannot exceed the SFU line. The 14.8× step between bf16 tensor cores and fp32 CUDA cores is the single largest cliff on the chip, and falling off it is entirely invisible in a FLOP count — the operations are identical, only the units executing them changed.

Two more effects sit on top. Clocks are not constant: sustained dense tensor-core work draws enough power that Hopper and Blackwell parts run below their boost clocks, so the real π during a long GEMM is meaningfully under the number on the slide. And a FLOP count carries a convention — a fused multiply-add is two operations here and throughout, which is near-universal, but vendor headline figures sometimes assume 2:4 structured sparsity and are then exactly twice the dense number.

None of this breaks T = W/π. It means π must be chosen honestly: the ceiling your instruction mix can actually reach, on the part you are actually running, at the clock it actually sustains.

Figure 02

Stage 01 / 07

Figure 02 — Which ceiling are you under? Compute ceilings are set by instruction mix; bandwidth slopes are set by access pattern. A gather-heavy fp32 kernel is bounded by two different degraded limits at once, and the roofline drawn against spec numbers would badly misdiagnose it.

Approximation · measure your own ceilings before trusting them
Every number in the stack above is a burst peak under ideal conditions, and β is no better — real coalesced streaming typically lands near 85% of published bandwidth. Comparing a measured point against spec roofs therefore systematically overstates your headroom, and the overstatement is largest exactly where you are already doing well. An empirical roofline, with both rates measured on the machine you are actually using in the datatype and access pattern you are actually using, is worth more than a precise model built on vendor maxima.

Memory time, and across which boundary

The data side is symmetric. Let Q be the bytes an operation must move and β the rate at which they can be moved. Then:

Memory timeT_memory = Q / β

This looks like the same equation with different letters, and it is — but Q hides a question that W does not, and the question is the most important one in this whole subject.

Trap · "bytes moved" is meaningless without a boundary
W is a property of the algorithm: a 4096³ matrix multiply is 137 GFLOP no matter what hardware runs it. Q is not. Bytes moved between what and what? The same kernel moves one amount across the HBM interface, a much larger amount across L2, a larger amount again between shared memory and the registers, and possibly zero across NVLink. These are different numbers describing different physical wires, and each has its own β. A statement about Q that does not name a boundary is not a claim that can be checked.

The boundaries that matter on a modern accelerator, from the arithmetic outward:

Where bytes can cross registers ↔ arithmetic units
shared memory / L1 ↔ registers
L2 ↔ shared memory
HBM ↔ L2 ← the default boundary throughout this article
NVLink ↔ peer GPU
network ↔ another node

Unless a section says otherwise, Q means traffic across the HBM boundary and β means HBM bandwidth. That is the right default because HBM is usually the tightest boundary that a kernel's design can actually influence. The full multi-level treatment comes later, once the single-boundary model is solid enough to generalise.

Counting Q honestly

Even with the boundary fixed, there are three defensible answers, and confusing them causes more argument than anything else here.

Compulsory traffic is what an oracle with unlimited cache would move: every input read exactly once, every output written exactly once. It is a lower bound on Q, and therefore gives an upper bound on intensity. Tiled traffic is what a specific blocking scheme moves, counting the re-reads that a finite cache forces. Measured traffic is what a hardware counter records, which includes over-fetch, replays and speculation, and is always the largest of the three.

All three are legitimate; they answer different questions. What is not legitimate is quoting one and reasoning with another. The convention here is that an unqualified Q is the compulsory count, and anything else is named.

Aside · a sanity check that always works
For any kernel, ask: what is the smallest number of bytes an oracle with unlimited cache would have to move to compute this output? That is Q_compulsory, it takes about thirty seconds to work out from the tensor shapes, and it bounds everything. If your measured traffic is 5× that, you have found a re-read problem before opening a profiler.

Nuance · "the arithmetic intensity" is not a number
It is a number per boundary, per accounting policy. The same GEMM has a compulsory HBM intensity in the thousands, a tiled HBM intensity around 64, a measured HBM intensity somewhere between, and a different number again at L2 — four correct answers to four different questions. Almost every disagreement about whether a kernel is memory-bound dissolves once both parties say which boundary and which policy they meant.

Arithmetic intensity

Arithmetic intensity is defined as

I = W / Q [FLOP · byte⁻¹]

and the single most important thing to understand about it is that it is not a property of the algorithm alone. It is a property of an algorithm and a chosen memory boundary and a chosen implementation. The same matrix multiply has one intensity measured at the HBM boundary, a different one at the L2 boundary, and a third at the register-file boundary. Papers that quote "the arithmetic intensity of attention" without saying which boundary they mean are quoting a number that cannot be checked.

The convention in this article, unless stated otherwise: Q is traffic crossing the HBM/DRAM boundary, because that is the boundary that dominates on real transformer workloads and the one profilers report most reliably.

The memory wall, quantified

Machine balance I_* = π/β is the hardware's side of the equation. It has been getting worse for a decade, and the reason is structural: arithmetic throughput scales with transistor count and with datatype narrowing, while off-chip bandwidth scales with pin count, signalling rate and stack height. The former has been much easier to grow.

Figure 03

Stage 01 / 09

Figure 03 — The widening gap. Dense tensor peak (lowest natively-supported floating format) and HBM bandwidth, indexed to P100 = 1. Compute has outrun bandwidth by roughly an order of magnitude across seven generations, and machine balance I_* rose with it. Vendor peak figures; B200 numbers are as published and not independently measured.

**Table 1 — Ridge points.** I_* = π/β for the HBM boundary, by accelerator and datatype. Note that narrowing the datatype raises I_* whenever peak scales faster than bytes shrink.
Accelerator	β (HBM)	π bf16	I* bf16	π fp8	I* fp8	π fp32	I* fp32
A100 SXM 80GB	2.04 TB/s	312 TF/s	153	—	—	19.5 TF/s	9.6
H100 SXM5	3.35 TB/s	989 TF/s	295	1979 TF/s	591	66.9 TF/s	20
H200 SXM	4.80 TB/s	989 TF/s	206	1979 TF/s	412	66.9 TF/s	14
B200 (spec)	8.00 TB/s	2250 TF/s	281	4500 TF/s	563	80 TF/s	10
MI300X	5.30 TB/s	1307 TF/s	247	2615 TF/s	493	163 TF/s	31
L40S	0.86 TB/s	362 TF/s	419	733 TF/s	848	91.6 TF/s	106
RTX 4090	1.01 TB/s	165 TF/s	164	330 TF/s	328	82.6 TF/s	82

Two readings worth internalising. First, the L40S column: a card with strong tensor cores and modest bandwidth has a brutal ridge of 419 FLOP/byte, which is why inference-serving economics on such parts depend so heavily on batching. Second, the fp32 column: on an H100, non-tensor fp32 code only needs 20 FLOP/byte to be compute-bound. Ridge points are not a property of the chip; they are a property of the chip and the instruction mix you choose.

Nuance · the ridge moves under you
I* is not a property of the roofline model, it is a property of a specific part in a specific datatype, and it has roughly doubled every hardware generation because π grows faster than β. An operator with I = 200 was compute-bound on an A100 (I* = 153) and is memory-bound on an H100 (I* = 295) without a single line of it changing. The same happens within one part when you narrow the datatype: fp8 doubles π and leaves β alone, so it doubles the ridge to 591 and pushes operators that were comfortably compute-bound back down onto the diagonal. Quantisation does not only make kernels faster, it moves the goalposts they are measured against.

Deriving the roofline

Let a kernel perform W floating-point operations and move Q bytes across some specified memory boundary. Two independent resources are consumed, so two independent lower bounds on execution time T hold simultaneously:

Two independent lower bounds T ≥ W / π (the arithmetic units cannot go faster than peak)
T ≥ Q / β (the memory system cannot go faster than peak)

If the two resources operate concurrently — loads issue while arithmetic retires — then the binding constraint is whichever takes longer, and the best achievable time is their maximum:

T_min = max( W/π , Q/β )

Turn this into a throughput. Attainable performance P is work divided by time:

P = W / T ≤ W / max(W/π, Q/β) = min( π , β · W/Q ) = min( π , β·I ) where I ≡ W/Q

That is the roofline. It took four lines. Everything else in this article is either an interpretation of those four lines or an exercise in computing W and Q honestly for a real kernel.

The overlap assumption is doing real work

The max is an idealisation. It assumes memory traffic and arithmetic overlap perfectly: while one warp waits on a load, another warp's MMA instructions issue. When overlap fails completely — a strictly phased kernel that loads, then computes, then stores — the times add instead:

No-overlap model T_serial = W/π + Q/β ⇒ P_serial = ( 1/π + 1/(β·I) )⁻¹

This is a harmonic mean rather than a minimum, and it produces a rounded roofline: no sharp corner, and a peak that is only asymptotically approached. Note the worst-case penalty is exactly 2× — at the ridge point, where W/π = Q/β, the serial model predicts half the overlapped throughput. Real kernels live between these two curves, and the gap between them is a decent proxy for "how well does this kernel hide latency".

Figure 04

Stage 01 / 07

Figure 04 — max versus sum. Left: the same work scheduled with and without overlap. Right: the roofline each schedule implies. The serialized model is exactly 2× below the overlapped model at the ridge and converges to it in both asymptotic limits.

Nuance · why the inequality is the interesting part
Because the bound is a maximum, improving the non-binding resource changes nothing. Halving the FLOP count of a memory-bound kernel buys you exactly zero. This is the single most common way engineering effort is wasted in ML systems work, and the roofline exists to prevent it.

Nuance · the whole model rests on one assumption
T ≥ max(W/π, Q/β) rather than W/π + Q/β holds only if the machine can move bytes and do arithmetic at the same time. Modern GPUs largely can — asynchronous copy engines, multi-stage software pipelines and a memory system that runs independently of the tensor cores exist precisely to make that true. But it is an assumption, not a theorem, and where it fails the real curve rounds off toward the sum, which sits exactly 2× below the roof at the ridge and approaches it far from the ridge. Everything downstream inherits this: the roofline is an upper bound that a machine with poor overlap will never reach, and a kernel sitting a consistent factor of two low may be failing to pipeline rather than failing to be efficient.

The roofline equation and its geometry

Plot P = min(π, β·I) with both axes on a logarithmic scale. Take logs of the bandwidth branch:

log P = log β + log I

A straight line of slope exactly 1, with vertical offset log β. The compute branch, log P = log π, is horizontal. The two meet at the ridge. This is why the roofline is always drawn log-log: it is the only coordinate system in which the memory bound is a straight line, and in which changing the hardware's bandwidth slides the diagonal vertically without changing its shape.

Three geometric facts fall out immediately, and they are worth stating explicitly because each corresponds to a real optimisation lever:

Reading the geometry 1. Horizontal motion (raising I) helps only left of the ridge.
2. Vertical distance to the roof = the headroom you can actually recover.
3. Fraction of peak: P/π = min(1, I / I_*) — a single scalar.

That third identity is the most useful one in practice. It says that a memory-bound kernel's achievable fraction of peak is exactly the ratio of its intensity to the machine balance, independent of everything else. An RMSNorm at I = 2 on an H100 can reach 2/295 ≈ 0.68% of peak bf16, and a perfect implementation of it will still reach 0.68% of peak. That is not a bug in the kernel. It is arithmetic.

Figure 05

Stage 01 / 10

Figure 05 — Building the roof, then furnishing it. Kernel placements assume bf16 storage on an H100 at the HBM boundary. The camera zooms into the memory-bound shelf, where almost every non-GEMM transformer operator lives, then out to the compute-bound plateau occupied only by large matrix multiplies and long-context attention.

Trap · "we cut the FLOPs" is not a performance claim
"We cut the FLOPs by 40%" is not a performance claim. If the kernel was memory-bound, cutting FLOPs moves its point left along the diagonal — the same wall time, a lower fraction of peak, and a worse-looking profile. Sparse and low-rank methods run into this constantly: the FLOP savings are real, the wall-clock savings are not, because the weights still have to be read.

Nuance · log–log flatters everything
Both axes are logarithmic, so vertical distance to the roof compresses badly. A point that looks like it is nearly touching the roof can easily be 2× below it, and the no-overlap machine of §5 sits exactly 2× under the ridge — a gap that is barely visible here and worth an entire engineering quarter. Read the ratio from the readouts, never from the picture.

Lab · the roofline explorer

Before deriving intensities for specific operators, get a feel for the plane. Pick hardware, a datatype, and a memory boundary; then slide the intensity or click a kernel preset. The verdict line tells you, for the current point, exactly how much more intensity you would need to reach the ridge.

Lab 1 · interactive

Accelerator

Datatype / pipe

Memory boundary

Arithmetic intensity 2.0

machine balance I*295

attainable P—

fraction of peak—

regime—

time for 1 TFLOP—

—

Lab 1. L2 and shared-memory bandwidths are approximate aggregate figures and vary strongly with access pattern; treat those two roofs as indicative, not authoritative. NVLink figures are per-GPU unidirectional.

Approximation · the on-chip roofs
The HBM and NVLink slopes in this lab come from published figures. The L2 and shared-memory slopes do not: they are aggregate estimates derived from issue width and clock, and both vary by a factor of two or more with access pattern, partition locality and hit rate. Use them to reason about which boundary is likely to bind and by roughly how much; do not quote them as bounds. Where a decision depends on the answer, measure the level with ncu and substitute your own number.

GEMM from first principles

Take C = AB with A ∈ ℝ^M×K, B ∈ ℝ^K×N, C ∈ ℝ^M×N, all stored at s bytes per element.

W = 2MNK
Q_compulsory = s(MK + KN + MN)
I_compulsory = 2MNK / [ s(MK + KN + MN) ]

For the square case M = N = K = n this collapses beautifully:

I_compulsory = 2n³ / (3sn²) = 2n / (3s)

which for bf16 (s = 2) is simply n/3. Set that equal to the H100 ridge of 295 and solve: n ≈ 885. So the compulsory bound says a square bf16 GEMM becomes compute-bound at roughly n = 885, and grows more compute-bound linearly after that. A 4096³ bf16 GEMM has a compulsory intensity of 1365 FLOP/byte — 4.6× past the ridge, with enormous margin.

Notice what the asymmetry does. If K is huge but M and N are small — a tall-skinny reduction — the MK and KN terms dominate and I → 2MN/(s(M+N)), a constant in K. If M = 1, I → 2N/(s(1+N)) ≈ 2/s. The intensity of a matrix multiply is governed by its smallest dimension, which is the entire reason autoregressive decoding is slow.

Figure 06

Stage 01 / 09

Figure 06 — Where the re-reads come from. Eight output tiles, each requiring a full row-panel of A and a full column-panel of B. Every panel is fetched once per tile that needs it. The counter tracks bytes actually pulled across the boundary against the compulsory minimum; the ratio is the tiling's re-read factor.

Approximation · compulsory traffic is a floor, so intensity is a ceiling
Q_compulsory assumes every byte is read exactly once, which requires a cache large enough to hold the whole working set. Real Q is therefore ≥ Q_compulsory, and real I is ≤ the compulsory I. Every compulsory intensity in this article is an optimistic bound, and when a kernel misses its compulsory roofline the first hypothesis should be that it is re-reading, not that the model is wrong.

GEMM II · tiles set the intensity

The compulsory bound is optimistic because it assumes an oracle that touches every byte once. Real GEMM kernels cannot: C does not fit on chip, so the output is decomposed into tiles of shape T_m × T_n, and each tile independently streams the corresponding row-panel of A and column-panel of B.

Count it. There are (M/T_m)(N/T_n) tiles. Each reads T_mK + KT_n elements and writes T_mT_n:

Q_tiled = s · (M/T_m)(N/T_n) · (T_mK + KT_n + T_mT_n)
= s · MNK · (1/T_n + 1/T_m) + s·MN

Divide into W = 2MNK and let K be large so the +sMN term vanishes:

The central result of GEMM blocking I_tiled = 2 / [ s · (1/T_m + 1/T_n) ] → T / s for square tiles T_m = T_n = T

Read that carefully, because it is the most important formula in this section and it surprises people: the intensity of a tiled GEMM does not depend on the matrix size at all. It depends on the tile size. Doubling n does nothing; doubling T doubles your intensity. This is why every high-performance GEMM library is fundamentally a tile-size search, and why 128×256×64 and 256×128×64 keep appearing in CUTLASS configuration tables.

And now the problem. With bf16 and a 128×128 tile, I_tiled = 64 FLOP/byte — well under the H100's ridge of 295. Even a very aggressive 256×256 tile gives 128. The tile that would be needed to reach 295 is T = 590, requiring an accumulator block of 590² fp32 values ≈ 1.4 MB per tile, which vastly exceeds the 228 KB of shared memory per SM. By the two-level model, a bf16 GEMM should be memory-bound on an H100. It is not.

What rescues it

Three things, in order of importance:

L2. The H100 has ~50 MB of L2. When tiles are scheduled in a swizzled or "threadblock rasterisation" order, the B-panel a tile needs was very likely just read by a neighbouring tile and is still resident. The re-reads happen — but at L2 bandwidth, not HBM bandwidth. HBM traffic drops back toward Q_compulsory, so the honest model has two boundaries, not one. That is the hierarchical roofline, and we build it properly later.

Split-K and stream-K scheduling. These change the loop structure to trade a small amount of extra output traffic for far better SM occupancy on shapes where the tile count is not a clean multiple of the SM count.

Deeper hierarchy inside the tile. The 128×128 threadblock tile is further blocked into warp tiles and register tiles; the register-level intensity is what actually feeds the tensor cores at 989 TFLOP/s. The HBM roofline never sees this.

Figure 07

Stage 01 / 08

Figure 07 — Two different lines, two different variables. Left: compulsory intensity grows with n. Right: tiled intensity grows with T and is flat in n. The shaded band is the region a two-level model says is unreachable and that L2 in fact makes reachable.

Nuance · T/s is an asymptote, not the answer
I → T/s drops the output term sMN, which is negligible only when K is large. For a skinny GEMM — K a few hundred, M and N large, which is exactly what a projection layer at short sequence length looks like — the output write is a substantial fraction of Q and the achieved intensity falls well below T/s. The tile size sets the ceiling; the aspect ratio decides how close you get to it. The SMEM wall at T ≲ 240 is likewise a consequence of fp32 accumulators occupying most of the 228 KB budget, so it moves if the accumulator precision moves.

Shapes matter

The GEMM analysis assumed a square matrix. Almost no matrix in a transformer is square, and the departures are not cosmetic — they are the difference between the same layer running at 70% of peak and at 0.3% of it.

Return to the compulsory intensity for a general M × K by K × N product:

Compulsory intensity, general shape I = 2MNK / [ s·(MK + KN + MN) ]

Nothing in that expression cares about the total size. It cares about the ratios. Divide through by MNK to see why:

I = (2/s) · 1 / ( 1/N + 1/M + 1/K )

The intensity is governed by a harmonic mean of the three dimensions, which means the smallest dimension dominates. Make one of M, N, K small and the whole expression collapses to roughly (2/s)·min(M, N, K), no matter how enormous the other two are. A 65536 × 65536 by 65536 × 8 GEMM does an enormous amount of arithmetic and is hopelessly memory-bound at I = 8, because every weight it reads is used exactly eight times before being discarded.

That single fact organises the entire performance story of a transformer, because the small dimension is almost always the same one: the token count.

**Table Table 2 — the same weight matrix, three regimes.** A d × d projection with d = 8192 in bf16, on an H100 whose ridge is 295. Only the token count changes.
regime	tokens M	shape	I = 2Md / [s(2M + d)]	vs ridge 295	bound
training / prefill	32768	fat	3641	12.3× over	compute
prefill, short prompt	512	fat-ish	455	1.5× over	compute, barely
decode, batch 64	64	skinny	63	4.7× under	memory
decode, batch 1	1	vector	1.0	295× under	memory

The weights are identical in all four rows and the arithmetic per weight byte is not, and it is the only thing that changed. This is why prefill and decode behave like different workloads running on different machines despite executing the same matrices, and why the whole discipline of inference serving is fundamentally an exercise in manufacturing token count.

The extreme case: decode is a GEMV

At batch 1, the projection degenerates to a matrix–vector product. Every weight is read from HBM, used for exactly one multiply-add, and discarded:

GEMV intensity I = 2d² / [ s·(d² + 2d) ] → 2/s = 1.0 FLOP/byte in bf16

Independent of model size, layer width, and every implementation choice available to you. The consequence is a step time that is a pure bandwidth quotient — read the whole model, emit one token:

t_step ≈ s·P / β ⇒ a 70B model in bf16 is 140 GB / 3.35 TB/s = 41.8 ms, about 24 tokens/s

Batching is the only lever that touches the numerator. With B sequences decoding together, each weight read serves B tokens, so I ≈ 2B/s and the ridge is reached at

B_* = s·I_* / 2 = 295 on an H100 in bf16

which is the reason continuous batching exists, and — since narrowing the datatype scales s and π together — the reason quantising to fp8 does not lower it.

Trap · B* is not a batch size you can actually run
B_* = 295 is the batch at which weight reads are amortised enough to reach the ridge. It says nothing about whether that batch fits or whether you can keep it full. Ragged sequence lengths leave a nominal batch half-empty; the KV cache usually exhausts memory long before B_*; and prefill work interleaved with decode steals the very bandwidth the amortisation depends on. Treat it as the target the scheduler is aiming at, not a number to configure.

Nuance · the KV cache changes this story again
The amortisation above covers the weights. It does not cover the KV cache, which is read in full every step and grows with both batch and context. For a 70B model at 8k context the cache is about 320 KB per token, so 2.6 GB per sequence — and past roughly 53 concurrent sequences the cache, not the weights, is the majority of all HBM traffic. Decode attention over that cache has I = 2/s regardless of how you schedule it, so the sub-linear region arrives far earlier than B_* alone suggests.

Figure 08

Stage 01 / 07

Figure 08 — Free until 295. The weight matrix is streamed once per step regardless of batch. Throughput therefore rises almost linearly until intensity reaches the ridge, after which the tensor cores become the constraint and per-token latency begins to grow.

Lab · the GEMM workbench

Set the three dimensions, the datatype, the tile shape and how much of the panel re-read traffic you believe L2 absorbs. The bench reports all three accountings of Q, places the resulting point on the roof and predicts a wall time. The L2 absorption slider is the honest knob: at 0% you get the pessimistic two-level model, at 100% you get the compulsory bound, and reality for a well-tuned CUTLASS kernel on large shapes tends to sit at 80–95%.

Lab 2 · interactive

M 4096

N 4096

K 4096

Datatype

Accelerator

Tile T 128

L2 absorption 85%

W (FLOPs)—

Q compulsory—

Q tiled—

Q modelled—

I modelled—

fraction of peak—

predicted time—

—

Lab 2. Predicted time is the roofline lower bound W/min(π, βI) and therefore optimistic: it ignores wave quantisation, launch overhead, clock throttling and epilogue costs. Treat the number as a floor, not a forecast.

Nuance · the absorption slider is a knob, not a prediction
The lab models L2's effect as a single fraction of the tiled re-read traffic that never reaches HBM, which is a deliberate simplification of something genuinely hard: real absorption depends on tile scheduling order, threadblock rasterisation, how much of the 50 MB the concurrent tiles collectively touch, and how the L2 partitions are addressed. Slide it to see the corridor the true answer lives in, then measure dram__bytes.sum to find out where in that corridor your kernel actually sits. The corridor is the useful output here; a single point inside it would be false precision.

Elementwise operations

Every operation that touches each element a constant number of times and does a constant amount of arithmetic with it has an intensity that is a pure constant — no dimension appears anywhere. There is nothing you can do to a residual add to make it compute-bound. Write the general form:

Elementwise intensity I = f / [ s·(r + w) ]
f = FLOPs per element, r = tensors read, w = tensors written, s = bytes per element

The denominator is at least 2s for any op that reads one tensor and writes one. So for bf16 the ceiling on any single-input single-output elementwise kernel is f/4, and f is a small integer. Concretely:

**Table 3 — The floor.** bf16 storage, fp32 accumulation, one-pass implementations with intermediates in registers. Fraction of peak is I/I_* against the H100 bf16 ridge of 295.
Operation	FLOP / elt	tensors r+w	bytes / elt	I	max % of bf16 peak
bf16 → fp32 cast	0	2	6	0.00	0.00%
residual add (x + y → z)	1	3	6	0.17	0.06%
Adam parameter update	12	7 tensors	30	0.40	0.14%
SwiGLU (elementwise part)	~6	3	6	1.00	0.34%
RMSNorm	~5	2	4	1.25	0.42%
softmax (online, one pass)	~6	2	4	1.50	0.51%
LayerNorm (mean + var)	~9	2	4	2.25	0.76%
GELU (tanh approximation)	~12	2	4	3.00	1.02%
fused norm+SwiGLU+residual	~12	2	4	3.00	1.02%

Nothing in that table reaches 1.1% of peak. That is not an implementation failure — it is the definition of the operations. The correct response is not to optimise these kernels but to stop launching so many of them, which is what the next section is about.

The useful corollary is that for any memory-bound op you can predict wall time to within a few percent using nothing but a division:

t ≈ s·(r + w)·n / β

Take activations of shape (batch 8, seq 8192, d 8192) in bf16: n = 5.4×10⁸ elements, 1.07 GB per tensor. A single RMSNorm reads and writes it: 2.15 GB at 3.35 TB/s = 0.64 ms. A transformer layer contains roughly a dozen such passes; eighty layers of them is over half a second of pure data movement per forward pass, before a single matrix multiply has been issued. Measured against the same layer's GEMM time, that is not a rounding error.

Figure 09

Stage 01 / 06

Figure 09 — Zooming out until the ridge appears. The first four stages use an axis that runs to 4 FLOP/byte. The fifth rescales continuously to 320 so the H100 ridge fits on the page; every bar collapses to a sliver. That collapse is the correct intuition for how far these operators sit from the roof.

Approximation · the FLOP counts in Table 3 are implementation-dependent
Whether GELU is the erf form or the tanh approximation, whether RMSNorm uses rsqrt or a square root and a divide, and whether the compiler contracts multiply-adds all shift f by a few operations. Individual rows in that table can be off by 30–50%. The conclusion is completely insensitive to it: every entry is below 3 FLOP/byte against a ridge of 295, so a 50% error in f moves nothing from one side of the ridge to the other. Precision in the numerator is not worth chasing when the verdict is decided by two orders of magnitude.

Fusion

Consider a fragment of a transformer block written the way a naive framework executes it: residual add, RMSNorm, elementwise scale, SiLU, elementwise multiply, dropout. Six kernels, six launches, and — this is the point — each one reads its input from HBM and writes its output back to HBM, because a kernel boundary is a memory boundary.

Let the tensor have n elements at s bytes. Unfused traffic is:

Q_unfused = Σ_i s·n·(r_i + w_i) ≈ 2s·n·(number of kernels) + extra for multi-input ops

Fuse them into a single kernel and every intermediate stays in registers. Traffic collapses to reading the genuine inputs and writing the genuine outputs:

Q_fused = s·n·(r_external + w_external)
I_fused = ( Σ_i f_i ) / [ s·(r_ext + w_ext) ] — numerator adds, denominator does not

Fusion is the cleanest example of intensity engineering in the whole stack: it leaves W exactly unchanged and divides Q by roughly the number of kernels fused. On the roofline the point translates horizontally, and because you were on the diagonal, performance rises by exactly the same factor.

What limits fusion

Three things, and they are all capacity or dependency constraints rather than bandwidth ones. Register and shared-memory pressure: a fused kernel must hold every live intermediate on chip, and spilling to local memory silently reintroduces the traffic you were trying to remove. Reductions: a normalisation needs a full row before it can scale anything, so fusing across a reduction requires either two passes or an online formulation. Training: the backward pass needs some intermediates, so fusing forward kernels forces a choice between storing them (traffic) and recomputing them (FLOPs) — and the roofline tells you which side of that trade you are on.

Figure 10

Stage 01 / 07

Figure 10 — Kernel boundaries are memory boundaries. Each fusion step removes one HBM round trip. FLOPs are constant across every stage; only the vertical arrows disappear. The intensity readout climbs by the same factor the traffic falls.

Trap · fusion can win on the roofline and lose on the clock
Fusing raises I exactly as the arithmetic says, but a fused kernel holds more live values, and register pressure costs occupancy. Below the concurrency Little's law demands (§17), the fused kernel cannot saturate β and gives back the bandwidth the fusion just saved — or worse, spills to local memory and adds traffic the model never counted. The roofline predicts the ceiling correctly and is silent about whether you can reach it. Fuse, then measure occupancy and spill counts before believing the prediction.

The big-vocabulary logits crisis

Language and vision-language models have drifted toward vocabularies of 100k–256k tokens, and this creates a roofline pathology at the very last layer that is worth working through in full because it is so often the largest single memory consumer in a training step.

Let N_t = B·T be the number of tokens in the microbatch, d the model width, V the vocabulary. The LM head is a GEMM of shape (N_t×d)(d×V):

W_head = 2·N_t·d·V
Q_head = s(N_td + dV + N_tV)

With N_t = 8192, d = 4096, V = 100,000 and bf16 this is 6.7 TFLOP against 2.5 GB — intensity 2660, comfortably compute-bound. So far so good. The problem is not the GEMM; it is the tensor it produces. The logits are N_t×V = 8.2×10⁸ values. In fp32, as cross-entropy conventionally wants them, that is 3.3 GB — and the loss needs at least the logits, the softmax output and the gradient, so the realistic footprint is 3–4 copies:

The pathology logits memory ≈ 4 · N_t · V bytes per copy
cross-entropy traffic ≈ 4 · 4 · N_t · V bytes (read logits, write softmax, read back, write grad)
I_CE ≈ 6 FLOPs / 16 bytes = 0.38 FLOP/byte

So the cross-entropy computation moves 13 GB for this microbatch — five times the traffic of the GEMM that produced the logits — at an intensity of 0.38, meaning roughly 0.13% of peak. It also, incidentally, consumes more activation memory than the entire rest of the network. On a 7B model with 100k vocabulary this single operator routinely dictates the maximum microbatch size.

The fix is a roofline fix

Chunked or fused cross-entropy processes the tokens in blocks of rows: for each chunk, compute its slice of logits, immediately reduce it to a loss contribution and a gradient with respect to the hidden state, and never materialise the full N_t×V tensor at all. The peak logits memory falls from 4N_tV to 4·C·V for chunk size C, and the entire cross-entropy traffic term disappears into on-chip work amortised inside the GEMM's epilogue.

Read the move on the roofline: the GEMM's W is unchanged, its Q loses the N_tV output term entirely, and a separate 0.38-intensity kernel vanishes. You do a small amount of recomputation in the backward pass to get there — the same trade FlashAttention makes, applied at the other end of the network.

Figure 11

Stage 01 / 06

Figure 11 — The last layer eats the budget. N_t = 8192, d = 4096. The vocabulary sweep runs from 32k to 256k. Chunked cross-entropy at C = 1024 cuts the peak logits memory by N_t/C and removes the low-intensity kernel from the traffic ledger entirely.

Trap · averaging intensity across a boundary hides the whole problem
The head GEMM sits at I ≈ 2660 and the cross-entropy at I ≈ 0.38. Report them together as one "loss layer" and you get an aggregate intensity in the hundreds that looks perfectly healthy and is completely uninformative — the compute-bound half drags the average across the ridge while the memory-bound half spends the time. Arithmetic intensity does not average meaningfully, because time does not: the slow operator dominates the sum while contributing least to the mean. Always decompose to the level at which one kernel has one bound.

Lab · the fusion bench

Toggle each boundary in the chain. A boundary that is on means the intermediate is written to HBM and read back by the next kernel; a boundary that is off means it stays in registers. Note that turning boundaries off never changes the FLOP count — only the traffic, and therefore only the time.

Lab 3 · interactive

Tokens (B×T) 65536

d_model 8192

Activation dtype

Accelerator

kernels launched6

HBM traffic—

FLOPs—

intensity—

time (roofline)—

speed-up vs all-split1.00×

—

Lab 3. Register pressure is modelled crudely: the bench warns when a fused group would need more than roughly 64 KB of live state per 256-thread block, which is where real fused kernels start spilling.

Nuance · what fusion can and cannot reach
Fusion removes the intermediate round trips between the operators you fuse. It cannot remove the traffic at the boundaries of the fused region — the input still has to be read once and the output written once — so the fused intensity converges to f_total/[s(r_ext + w_ext)] and no further. That ceiling is usually still deep in memory-bound territory: fusing every elementwise operator in a transformer block gets you to about 3 FLOP/byte against a ridge of 295. Fusion is worth doing because it removes a large constant factor, not because it changes which side of the ridge you are on. Expecting it to make elementwise work compute-bound is expecting the wrong thing from it.

Attention from first principles

Attention is where all of this stops being an exercise. It is the operation whose performance story cannot be told in FLOPs at all, and the one where a data-movement analysis produced an algorithm that changed how every model is trained.

Trace it honestly, one tensor at a time, for a single head with sequence length N and head dimension d. Start from the projections: Q, K and V are each an N × d output of an N × d by d × d GEMM. Those are ordinary GEMMs and behave exactly as the previous sections describe — fat when N is large, degenerate when N is 1.

Then the part that is not a GEMM:

The core of attention S = QKᵀ      (N × d) · (d × N) → N × N
P = softmax(S)      N × N → N × N
O = PV      (N × N) · (N × d) → N × d

Ask the question that matters: which of these intermediates gets materialised, and how large is it? Q, K, V and O are all N × d. S and P are N × N. At N = 4096 and d = 128, the inputs are half a megabyte each and the intermediates are sixteen megabytes each — a factor of thirty-two, growing linearly with sequence length and entirely absent from any parameter count.

Now count the traffic of the naive implementation, in which S and P are real tensors in HBM. S is written and read back for the softmax; P is written and read back for the second GEMM. That is four N² round trips against 4N²d FLOPs:

W = 4N²d Q = 4s(Nd + N²) ⇒ I = Nd / [s(N + d)] → d/s = 64 as N grows

An intensity that converges to a constant set by the head dimension, against a ridge of 295. Naive attention is capped at 21.7% of peak and gets no better with longer sequences — the arithmetic grows as N², but so does the traffic, and the ratio flattens. Growing the problem does not help, because the problem was never the arithmetic.

Nuance · what the N² accounting does and does not assume
The count of four N² round trips assumes S and P are separately materialised, which is what an unfused implementation does. Fusing the softmax into the score kernel removes one pair and improves the constant without touching the asymptote — I still tends to d/s, because the N² traffic scales with N² exactly as the arithmetic does. Note also that d/s is a limit as N grows: at short sequence lengths the Nd terms are not negligible and the true intensity is higher than 64. The plateau is what long context converges to, not what every naive kernel measures.

FlashAttention is an IO algorithm

FlashAttention tiles the computation so that S and P exist only in on-chip memory. Split Q into blocks of B_r rows and K, V into blocks of B_c rows. For each Q block, stream all K/V blocks through shared memory, maintaining a running maximum and a running normaliser so that softmax can be computed online without ever seeing the full row at once. The output block is rescaled as each new K/V block arrives.

Nothing about the arithmetic changes. W is still 4N²d — in fact the backward pass does more arithmetic, because S is recomputed rather than stored. What changes is Q.

The compulsory bound

If K and V could be held on chip for the whole kernel, each of Q, K, V would cross the boundary exactly once and O would be written once:

Q_compulsory = 4sNd ⇒ I = 4N²d / (4sNd) = N / s

Linear in sequence length, independent of head dimension. At N = 2048 and bf16 that is I = 1024, more than three times the H100 ridge. Long context is, in principle, deeply compute-bound. This single fact is why context lengths grew as fast as they did once the IO-aware kernel existed.

The tiled bound

They cannot all be held on chip. FlashAttention-2 loops over Q blocks on the outside and re-reads the entire K and V for each one:

Q_tiled = s·Nd (Q) + s·(N/B_r)·2Nd (K,V re-read per Q block) + s·Nd (O)
I_tiled = 4N²d / [ s(2Nd + 2N²d/B_r) ] → 2B_r / s for N ≫ B_r

With B_r = 128 and bf16 that is 128 FLOP/byte — the exact same structural answer as GEMM tiling, and for exactly the same reason. Once again the truth is in between: K/V blocks re-read by consecutive Q blocks are usually still in L2, so measured HBM traffic sits far closer to the compulsory bound than to the tiled bound.

Nuance · what FlashAttention actually changed
FlashAttention is not a FLOP optimisation, an approximation, or a numerics trick. It is a change to Q that leaves W alone, which on the roofline is a purely horizontal move. Recomputation in the backward pass is the same move in a more aggressive form: deliberately increase W to decrease Q, because the point was to the left of the ridge and the exchange rate was favourable.

Figure 12

Stage 01 / 10

Figure 12 — Same FLOPs, different Q. Top: the N×N score matrix round-tripping through HBM four times. Bottom: the tiled kernel keeping every intermediate on chip and rescaling the accumulator with an online softmax. The counters use N = 4096, d = 128, bf16, B_r = B_c = 128.

Figure 13

Stage 01 / 06

Figure 13 — Three curves, one operator. Naive attention plateaus at d/s. The tiled Flash kernel plateaus at 2B_r/s. Only the compulsory accounting grows without bound. Real kernels track the compulsory line closely at long context because L2 captures the K/V re-reads.

Nuance · causal masking is not an optimisation here
Restricting attention to the lower triangle halves W and halves Q, so I is unchanged and the point does not move on the roofline at all. The kernel gets faster because there is less of everything, not because it got more efficient — its position relative to the roof is identical. This is the cleanest example of why "we cut the FLOPs" and "we improved the intensity" are different claims: only the second one moves you up.

Attention III · KV cache, GQA, MLA

Decoding inverts everything. The query is a single vector; K and V are the accumulated cache of length N. Per head:

W = 2Nd (qKᵀ) + 2Nd (pV) = 4Nd Q = 2s·Nd (read K and V)
I_decode = 4Nd / (2sNd) = 2/s = 1 FLOP/byte

Independent of N. Independent of d. And — crucially — independent of batch size, because every sequence owns a private cache. Batching amortises weights; it does not amortise KV. Decode attention is permanently pinned at the far-left edge of the roofline.

Grouped-query attention is the one structural lever. If h query heads share h_kv key/value heads with group size g = h/h_kv, the same K/V bytes feed g times as much arithmetic:

I_GQA = 2g / s (MQA is the extreme g = h)

Llama-3-70B uses h = 64, h_kv = 8, so g = 8 and I = 8 — an 8× improvement that is still 37× short of the ridge, but that also shrinks the cache by 8×, which is where the real win lives. Multi-head latent attention (MLA) goes further by caching a low-rank latent and reconstructing K/V on the fly: it deliberately spends FLOPs to buy back bytes, which is a rightward move on exactly the same plot.

The crossover that ends batching

Per sequence, the cache costs

KV bytes = 2 · s · N · d_head · h_kv · L

For Llama-3-70B (L = 80, h_kv = 8, d_head = 128, bf16) that is 320 KB per token, so a 8192-token context is 2.6 GB per sequence. Total decode traffic per step is weights plus caches:

Q_step = s·P_params + B · (KV bytes) ⇒ B_cross = s·P_params / (KV bytes)

Here B_cross = 140 GB / 2.6 GB ≈ 53. Past a batch of ~53 sequences at 8k context, more than half of all HBM traffic is KV cache, and the sub-linear region of the throughput curve arrives far earlier than the naive B_* = 295 suggested. This is the quantitative reason paged caches, cache quantisation, GQA and MLA are all load-bearing rather than nice-to-have.

Figure 14

Stage 01 / 06

Figure 14 — Weights are shared, caches are not. Per-step HBM traffic for a 70B model at 8k context, bf16. The weight term is constant in B; the cache term is linear. Reducing h_kv moves the crossover right and is the reason every recent model ships with grouped or latent attention.

Nuance · this one cannot be fixed in the kernel
Decode attention has I = 2/s independent of N, of B, and of every scheduling decision available to you. No tiling, fusion, launch configuration or memory layout changes it, because the arithmetic per KV byte is fixed by the algorithm. Every real remedy therefore changes the model rather than the code: fewer KV heads (GQA, MQA), a compressed cache (MLA), fewer cached bytes per entry (quantisation), or fewer entries (windowing, eviction). That is an unusual and important situation — the roofline is telling you the bottleneck is not an implementation problem.

Lab · the transformer layer budget

This is the synthesis. Every operator in one transformer layer, with its own W, its own Q, its own intensity, and its own roofline-predicted time. The stacked bar shows where the time actually goes; the MFU readout tells you what fraction of the machine you are using. Switch between prefill, a training step and decode, and watch the entire shape of the budget invert.

Lab 4 · interactive

Regime

d_model 8192

Sequence T 8192

Batch B 4

KV heads

Datatype

Accelerator

layer FLOPs—

layer HBM traffic—

layer time—

MFU (this layer)—

time memory-bound—

—

Lab 4. Head dimension is fixed at 128 and d_ff at 3.5·d. Training multiplies FLOPs by 3 and traffic by 2 as a first-order stand-in for the backward pass. Times are roofline lower bounds and sum as if kernels ran back to back with no overlap between them.

Aside · three things worth trying in Lab 4
Switch to decode and watch every row turn red — that is the entire inference-serving problem in one screen. Turn off FlashAttention at T = 32768 and watch the attention row grow by two orders of magnitude. Drop d_model to 1024 in prefill and watch the GEMMs fall off the compute roof, which is why small models are so much harder to run efficiently than large ones.

Approximation · this budget is a bound on a bound
Three separate optimisms compound here. Each operator's time is its own roofline lower bound, so each row is already optimistic. The rows are then summed as if kernels ran strictly back to back, which ignores both the launch gaps that inflate the real total and the overlap that deflates it. And the training mode multiplies forward FLOPs by 3 and traffic by 2 as a stand-in for a backward pass whose real traffic depends on what you recompute versus what you stash. The relative shape — which operators dominate, which side of the ridge they sit on, how that inverts between prefill and decode — is robust to all three. The absolute milliseconds are not; expect real kernels to land 1.3–2× above these figures.

Optimizers as memory workloads

One more workload deserves a mention here, not because it is subtle but because it is routinely left out of performance accounting altogether: the parameter update.

Count Adam with fp32 master weights, per parameter. It reads the fp32 master weight, the gradient, and the two moment estimates, then writes back the two moments, the master weight, and the bf16 parameter the forward pass will use:

Adam traffic per parameter read master 4 B + grad 4 B + m 4 B + v 4 B = 16 B
write m 4 B + v 4 B + master 4 B + param 2 B = 14 B
30 bytes moved to perform about 12 FLOPs

I_Adam = 12 / 30 = 0.40 FLOP/byte ⇒ 0.14% of H100 bf16 peak

Nothing in training is more purely memory-bound. For a 7B model on one H100 that is 210 GB at 3.35 TB/s = 62 ms every step, independent of batch size and sequence length. On a 400 ms step, Adam alone is 15% of the wall clock while contributing a rounding error to the FLOP count — which is exactly how it disappears from MFU-based analyses.

The lever is visible directly in the byte count, and it is the state, not the arithmetic. Plain SGD moves 14 B/param, SGD with momentum 22 B, Adafactor roughly 22 B with its factored second moment, 8-bit Adam about 16 B. Each is a proportional reduction in optimizer step time, because every one of them sits so far below the ridge that time is exactly proportional to bytes.

The interesting exception is Muon, whose Newton–Schulz orthogonalisation runs five iterations of a degree-5 polynomial in the momentum matrix, three n × n matmuls apiece. For an n × n weight that is ≈30n³ FLOPs against ≈4n² bytes of state, giving I ≈ 7.5n/s — thousands of FLOP/byte. Muon is the only widely used optimizer whose step is compute-bound, which is precisely why it can be so much more expensive per parameter than Adam while costing little additional wall time: it spends a resource that was sitting idle.

The optimization theory behind these updates — why the moments exist, what Newton–Schulz is orthogonalising, how the sharding works — is the subject of The anatomy of a step. What matters here is only the byte count, and the byte count says that choosing an optimizer is partly a bandwidth decision.

Nuance · sharding relocates the traffic, it does not delete it
Under ZeRO-2 or FSDP each rank reads 30·P/world optimizer bytes instead of 30·P, and the HBM cost falls proportionally. The bytes have not disappeared; they have moved to a different boundary, where they reappear as reduce-scatter and all-gather traffic governed by a roofline whose β is the interconnect. Whether that is a win depends on whether β_net·world exceeds β_HBM — which is why the same sharding decision is obviously right inside a node and often wrong across a slow one, and why it needs the next article rather than this one.

Figure 15

Stage 01 / 07

Figure 15 — Seven streams for twelve FLOPs. The Adam update reads four tensors and writes three, doing barely a dozen operations with each parameter. The lower panel converts bytes into wall time for a 7B model and shows what sharding across ranks does to it.

The third regime: latency

The roofline has a blind spot that trips up small kernels constantly. Both bounds assume the relevant resource can be saturated, and saturating a memory system requires enough independent work in flight. Little's law makes this precise. For any queueing system in steady state:

Little's law, applied to memory concurrency = throughput × latency ⇒ bytes in flight = β × t_latency

On an H100 with β = 3.35 TB/s and an HBM load-to-use latency of roughly 600 ns, the machine must keep

3.35 × 10¹² B/s × 600 × 10⁻⁹ s ≈ 2.0 MB in flight

continuously in order to reach peak bandwidth. Spread over 132 SMs that is ~15 KB per SM, or roughly 120 outstanding 128-byte requests per SM. A kernel that cannot generate that much memory-level parallelism will not reach β regardless of its intensity, and the achieved bandwidth is instead

β_eff = min( β , C / t_latency ) where C is bytes actually in flight

This creates a third regime, latency-bound, that sits underneath the diagonal rather than on it. It is why a kernel launched with 4 thread blocks runs at a fraction of peak bandwidth even though it is nominally "memory-bound and therefore optimal", and it is the mechanism behind most of the mystery in small-batch inference.

Two ways to generate the concurrency, and they are interchangeable. Occupancy: more warps per SM, each with one load outstanding. Instruction-level parallelism: fewer warps, each with several independent loads outstanding — the classic result that a well-unrolled kernel at 25% occupancy can match a naive one at 100%. Little's law does not care which you use; it only counts bytes in flight.

Figure 16

Stage 01 / 06

Figure 16 — Filling the pipe. Left: outstanding requests accumulating toward the 2 MB the machine needs in flight. Right: achieved bandwidth against concurrency, and the resulting latency-bound shelf below the roofline diagonal.

Trap · a point below the diagonal is not evidence of a bad kernel
The roofline's two bounds are both saturation bounds, so the model has no vocabulary for "not enough work in flight" and simply draws you underneath everything. Before concluding that a kernel needs algorithmic work, check the three things the model cannot see: is the grid large enough to fill every SM, are there enough independent loads per thread to cover 600 ns, and is the kernel long enough that the launch is not a material fraction of it? A great many kernels that look like roofline failures are Little's-law failures, and the fixes are entirely different.

The hierarchical roofline

Everything so far used a single memory boundary. Real machines have five or six, each with its own bandwidth, and a kernel has a different arithmetic intensity at every one of them. The honest model is therefore not one roof but a family:

Multi-level bound P ≤ min( π , β_smem·I_smem , β_L2·I_L2 , β_HBM·I_HBM , β_net·I_net )

Each level has its own ridge π/β_ℓ, and the level that binds is whichever term in that minimum is smallest. This resolves the GEMM paradox from §7 cleanly. A 4096³ bf16 GEMM with 256×128 threadblock tiles has I_L2 ≈ 85 against an L2 ridge of ~110, so L2 is the tightest on-chip boundary and caps it at roughly 77% of peak — which is about where large GEMMs actually measure. At the shared-memory level its 64×64 warp tiles give I_smem ≈ 32 against a ridge of ~33, so SMEM is marginal but adequate. At HBM its intensity is near the compulsory 1365 against a ridge of 295, so HBM is not remotely the constraint. Three boundaries, three different verdicts, all of them correct.

**Table 4 — H100 boundaries.** Ridge is π_bf16/β. Register-file, shared-memory and L2 figures are order-of-magnitude aggregates derived from issue width and clock, not vendor-published bandwidths; treat them as approximate.
Boundary	β (approx.)	ridge I*	capacity	what it limits
register file	~300 TB/s	~3	64 K × 32-bit / SM	MMA operand supply
shared memory	~30 TB/s	~33	228 KB / SM	warp-tile reuse
L2 cache	~9 TB/s	~110	50 MB	threadblock-tile reuse
HBM3	3.35 TB/s	295	80 GB	weights, activations, KV
NVLink 4 (per GPU)	450 GB/s	2198	—	TP, intra-node FSDP
PCIe Gen5 ×16	64 GB/s	15459	—	host offload
EFA (p5, per GPU)	50 GB/s	19788	—	inter-node collectives

The ridge column is the useful one. It says that to be compute-bound across NVLink you need 2200 FLOPs per byte crossing the link, and across a 400 Gb/s network interface you need nearly 20,000. Those numbers govern how parallelism strategies behave, which is §20.

Figure 17

Stage 01 / 08

Figure 17 — One kernel, five verdicts. A 4096³ bf16 GEMM plotted at each boundary. Filled markers are boundaries its data actually crosses; hollow markers show where it would sit if every tile had to travel over the interconnect, which is why GEMMs are not split across devices this way. Only on-device boundaries are eligible to bind, and here L2 is the one that does.

Nuance · "memory-bound" is not a property of a kernel
It is a property of a kernel at a boundary. The 4096³ GEMM in Figure 17 is simultaneously compute-bound at HBM, L2-bound on chip, and marginal at shared memory — all three at once, all three true. The single most common analysis error in this area is carrying a verdict established at one level across to another: concluding from "the GEMM is compute-bound" that its L2 behaviour does not matter, or from "attention is memory-bound" that a bigger tile cannot help. Name the boundary every time you name the bound.

Performance diagnosis

The model earns its keep as a diagnostic. Given a kernel that is slower than you want, it turns an open-ended question into a short sequence of closed ones, each with a different answer and a different fix.

You need three numbers. W, counted analytically from the shapes as we have done throughout — not read from a profiler, which frequently miscounts tensor-core work. Q, from hardware byte counters at the boundary you care about. t, from a timing loop with warmup and enough iterations that launch overhead is separable. From those, I = W/Q and P = W/t, and the point lands somewhere relative to the roof.

The three positions On the diagonal — bandwidth-limited, running well. The fix is to reduce Q.
On the ceiling — compute-limited, running well. The fix is to reduce W or narrow the datatype.
Below both — something the model cannot see. Do not touch the algorithm yet.

The third case is the common one and the one worth spending time on, because "below both" is not a diagnosis, it is a prompt for a second round of questions.

The questions, in order

Are the tensor cores being used at all? A kernel running on fp32 CUDA cores is measured against a ceiling 14.8× lower than the headline. Before concluding a kernel is inefficient, establish which ceiling it is under. A great many "we are at 5% of peak" reports are 70% of the correct, much lower, ceiling.

Is HBM actually saturated? Compare measured bandwidth against a STREAM-like measurement on the same machine, not against the spec sheet — real coalesced streaming typically lands near 85% of published β, and uncoalesced access can lose another 8× to sector over-fetch. A kernel at 40% of spec bandwidth may be at 90% of what its access pattern permits, in which case the fix is the access pattern rather than the kernel.

Is there enough parallel work in flight? Both roofline bounds are saturation bounds, so the model simply has no vocabulary for "not enough work". This is the subject of the latency section and the most frequent reason a well-written kernel sits below the diagonal.

Are the kernels too small? A launch costs a few microseconds. Any kernel with a predicted time under ~20 µs is partly measuring its own launch, and a chain of a hundred tiny elementwise kernels is dominated by it. Check the kernel duration histogram before the kernel internals.

Is reuse happening where you assumed? The gap between analytic Q and measured Q is the most informative single number in this exercise. Measured Q is always ≥ algorithmic Q; when it is far larger, the kernel is re-reading data you assumed was cached, and the size of the gap is the size of the prize.

Is the work quantised badly? Thread blocks are dispatched in waves across SMs. With 132 SMs and 140 tiles you run two waves, the second 6% occupied, and the kernel takes twice as long as its FLOP count implies. The roofline is a function of totals and cannot see this at all.

Answering these in order is deliberate: each one is cheaper to check than the one after it, and each can fully explain a result that would otherwise send you rewriting an algorithm that was never the problem. The tooling that produces the numbers — counter names, profiler configuration, timeline analysis — is worth a treatment of its own and gets one later in this series.

Nuance · measured intensity is a lower bound on algorithmic intensity
Hardware byte counters record what the memory system actually moved, which includes sector over-fetch on uncoalesced access, replayed instructions, speculative prefetch and read-for-ownership on partial writes. Measured Q is therefore ≥ the Q your algorithm requires, and measured I ≤ the I you derived on paper. That gap is not noise, it is the signal, and its size is the size of the prize.

Figure 18

Stage 01 / 06

Figure 18 — From counters to a verdict. Three measurements produce one point, and the point produces one of three actions. The dashed vertical is the distance to the roof — the only headroom that actually exists.

What the model still cannot see

Being explicit about the failure modes is what separates using this model from believing it. It cannot see wave quantisation or launch overhead, as above. It cannot see that caches are modelled rather than simulated — Q depends on size, layout and scheduling in ways no formula predicts a priori, which is why measured intensity should always beat predicted intensity when they disagree. It has no notion of dependencies: a real step is a graph of kernels with a critical path and CPU-side gaps, so per-kernel bounds sum to a bound on the step, but a loose one. And it ignores asymmetries — reads and writes do not cost the same, read-for-ownership turns partial writes into read-modify-writes, and chiplet parts add non-uniformity that a single β cannot express.

None of this makes the model less useful. It makes it a first question rather than a last answer.

The next wall

Every boundary in this article has been inside one GPU. The model never assumed that: β was only ever "the rate at which bytes cross the boundary in question", and nothing about the derivation cared whether that boundary was a memory controller or a cable.

So point it outward. Replace β with the bandwidth of an interconnect and Q with the bytes a rank sends or receives, and the same inequality produces a ridge for distributed execution:

The same model, one level out I_net = FLOPs per byte crossing the link ⇒ communication hides behind compute ⟺ I_net ≥ π / β_net

The ridge points that come out are of a completely different order from anything we have seen:

HBM3, 3.35 TB/s   → I* = 295
NVLink 4, 450 GB/s per GPU   → I* = 2,198
PCIe Gen5 ×16, 64 GB/s   → I* = 15,459
400 Gb/s fabric, 50 GB/s per GPU   → I* = 19,788

Those numbers are the whole reason distributed training looks the way it does. Nearly twenty thousand FLOPs per byte is a bar that almost nothing clears by accident, which is why some parallelism strategies are confined to a single node while others cross datacentres, and why the arrangement of GPUs into nodes is a performance decision rather than a packaging one. Working that out properly needs its own treatment — collective algorithms, overlap, the arithmetic of each parallelism strategy — and that is the next article.

Nuance · the ridge is the only thing that changed
It is worth being clear about what this section is and is not claiming. The inequality is identical, the derivation is identical, and arithmetic intensity means exactly what it has meant throughout — FLOPs per byte crossing a boundary. All that moved is which boundary, and therefore how large I* is. That is why the network case belongs in a separate article rather than a longer version of this one: the model needs no extension, only new numbers and the collective algorithms that produce them.

Everything up to here assumed the bytes you need already live on this GPU. The remaining question is what happens when they do not.

References

Grouped by topic. Every entry links to its primary source. Vendor peak figures quoted in this article come from the architecture whitepapers listed under Hardware; L2, shared-memory and register-file bandwidths are approximations derived from issue width and clock, not published numbers, and are labelled as such wherever they appear.

The model itself

Roofline: An Insightful Visual Performance Model for Multicore Architectures · Williams, Waterman & Patterson, CACM 2009 — the original paper; introduces ceilings, walls and the ridge point.
Cache-Aware Roofline Model: Upgrading the Loft · Ilic, Pratas & Sousa, IEEE CAL 2014 — reformulates intensity from the core's perspective so cache levels appear as distinct roofs.
Hierarchical Roofline Analysis for GPUs · Yang, Kurth & Williams, CCPE 2019 — the multi-level GPU methodology used in §18.
Hierarchical Roofline Analysis: How to Collect Data on Intel CPUs and NVIDIA GPUs · Yang et al., 2020 — the practical counter recipes behind §19.
Hierarchical Roofline Performance Analysis for Deep Learning Applications · Yang, Wang, Farrell & Williams, 2020 — the model applied specifically to DL kernels.
NERSC Roofline documentation · a maintained, practical guide to collecting roofline data.

Memory systems, latency and measurement

STREAM: Sustainable Memory Bandwidth in High Performance Computers · McCalpin — the origin of "machine balance" as a measured quantity.
Understanding Latency Hiding on GPUs · Volkov, UC Berkeley 2016 — the definitive treatment of Little's law on GPUs and why occupancy and ILP are interchangeable.
Dissecting the NVIDIA Volta GPU Architecture via Microbenchmarking · Jia et al., 2018 — measured latencies and bandwidths at every level of the hierarchy.
Nsight Compute Profiling Guide · NVIDIA — counter definitions including dram__bytes.sum and the roofline section files.

GEMM, tiling and tensor cores

CUTLASS · NVIDIA — the tile-shape taxonomy that makes I = T/s concrete.
Matrix Multiplication Background User's Guide · NVIDIA — tile quantisation, wave quantisation and the arithmetic-intensity discussion for GEMM.
Stream-K: Work-Centric Parallel Decomposition for Dense Matrix-Matrix Multiplication · Osama et al., 2023 — scheduling that removes the wave-quantisation tail.

Attention and IO-awareness

FlashAttention: Fast and Memory-Efficient Exact Attention with IO-Awareness · Dao et al., 2022 — the canonical statement that attention is a data-movement problem.
FlashAttention-2 · Dao, 2023 — the Q-outer / KV-inner loop order analysed in §14.
FlashAttention-3 · Shah et al., 2024 — asynchrony and fp8 on Hopper.
Self-attention Does Not Need O(n²) Memory · Rabe & Staats, 2021 — the online-softmax formulation.
Data Movement Is All You Need · Ivanov et al., 2020 — a full data-movement analysis of transformer training.
Fast Transformer Decoding: One Write-Head is All You Need · Shazeer, 2019 — MQA, and the clearest early statement of the decode bandwidth wall.
GQA: Training Generalized Multi-Query Transformer Models · Ainslie et al., 2023.
DeepSeek-V2 · DeepSeek-AI, 2024 — multi-head latent attention as a bytes-for-FLOPs trade.
Efficient Memory Management for LLM Serving with PagedAttention · Kwon et al., 2023.

Training systems

Distributed execution is the subject of the next article in this series; these are the two this one actually leans on.

Reducing Activation Recomputation in Large Transformer Models · Korthikanti et al., 2022 — selective recompute as an explicit trade of W for Q.
Training Deep Nets with Sublinear Memory Cost · Chen et al., 2016 — gradient checkpointing, the original version of that trade.

Fusion, kernels and the loss layer

Liger Kernel: Efficient Triton Kernels for LLM Training · Hsu et al., 2024 — fused RMSNorm, RoPE, SwiGLU and fused-linear cross-entropy.
Cut Your Losses in Large-Vocabulary Language Models · Wijmans et al., 2024 — the cross-entropy analysis behind §11.
Triton · OpenAI — the language most fused kernels in this space are now written in.

Optimizers and precision

Adam: A Method for Stochastic Optimization · Kingma & Ba, 2015.
Adafactor: Adaptive Learning Rates with Sublinear Memory Cost · Shazeer & Stern, 2018.
8-bit Optimizers via Block-wise Quantization · Dettmers et al., 2021 — state bytes as the direct lever on step time.
Muon: An Optimizer for Hidden Layers in Neural Networks · Jordan et al., 2024 — the primary write-up of the Newton–Schulz update.
Muon is Scalable for LLM Training · Liu et al., 2025 — distributed Muon at scale.
Mixed Precision Training · Micikevicius et al., 2018 — the origin of the master-weight byte accounting in §16.
FP8 Formats for Deep Learning · Micikevicius et al., 2022.

Hardware

NVIDIA H100 Tensor Core GPU Architecture · NVIDIA whitepaper — source for the H100 peak and HBM3 figures used throughout.
NVIDIA Blackwell Architecture · NVIDIA — B200 figures, quoted as published.
AMD Instinct MI300X · AMD — CDNA 3 peak and HBM3 bandwidth.
Elastic Fabric Adapter · AWS — per-GPU network bandwidth for the p5 family, quoted in §20.