Roofline model

Deep dive

How it works

In: the kernel's flop count and the bytes it moves, plus the machine's peak flop/s and measured bandwidth · Out: the attainable flop/s bound and which resource limits it · State: none · Cost: O(1) arithmetic per kernel.

1. Peak compute = cores × GHz × FMA units × SIMD lanes × 2 (an FMA is 2 flops). For example, 16 cores × 3.5 GHz × 2 FMA × 8 f32 lanes (AVX2) × 2 = 1792 GFLOP/s.
2. Bandwidth: measure it (STREAM triad), don't take the datasheet. Dual-channel DDR5-4800 is 2 × 38.4 = 76.8 GB/s nominal. Triad usually reaches order of 70-90 % of nominal.
3. Ridge = 1792 / 76.8 ≈ 23 flop/B.
4. AI of the kernel = flops / bytes that cross the memory level you care about (DRAM by default). Count the bytes the hardware moves, not the bytes the code names. A store to a line not in cache first reads it (write-allocate), unless you use non-temporal stores.
5. Bound = min(peak, AI × BW). Compare with the measured speed. Far below the bound means a latency or parallelism problem. At the bound, change the algorithm.

Worked example (machine above):

Tiled matmul: for each tile product, load one b×b tile of A and one of B (2 × 4b² bytes) and do 2b³ flops, so AI = b/4. With b = 128 the three tiles take 3 × 64 KiB = 192 KiB. That fits in L2, and AI = 32 exceeds the ridge. See tiled matmul.

LLM decode: batch-1 matrix-vector uses each BF16 weight (2 B) for 2 flops, so AI ≈ 1. A 7B-parameter model (14 GB) on 3.35 TB/s therefore tops out near 3350 / 14 ≈ 240 tokens/s, whatever the tensor-core peak. Batching B sequences multiplies AI by about B. This is why continuous batching pays.

Code

Measures bandwidth with a multithreaded STREAM triad and prints the bound for a few kernels.

use std::hint::black_box;
use std::time::Instant;

/// Measured DRAM bandwidth: STREAM triad a[i] = b[i] + s*c[i] on all cores.
/// Bytes per element = 3 x 8 (reads b, c; writes a) + 8 for the write-allocate read of a.
fn triad_gbs(n: usize, threads: usize, reps: usize) -> f64 {
    let mut a = vec![0.0f64; n];
    let b = vec![1.0f64; n];
    let c = vec![2.0f64; n];
    let chunk = n.div_ceil(threads);
    let mut best = f64::INFINITY;
    for _ in 0..reps {
        let t = Instant::now();
        std::thread::scope(|sc| {
            for ((a, b), c) in a.chunks_mut(chunk).zip(b.chunks(chunk)).zip(c.chunks(chunk)) {
                sc.spawn(move || {
                    for i in 0..a.len() { a[i] = b[i] + 3.0 * c[i]; }
                });
            }
        });
        best = best.min(t.elapsed().as_secs_f64());
    }
    black_box(&a);
    (n * 32) as f64 / best / 1e9
}

/// Roofline: attainable GFLOP/s = min(peak, AI x bandwidth).
fn attainable(peak_gflops: f64, bw_gbs: f64, ai: f64) -> f64 { peak_gflops.min(ai * bw_gbs) }

fn main() {
    let threads = std::thread::available_parallelism().map_or(1, |n| n.get());
    let bw = triad_gbs(1 << 25, threads, 5); // 3 arrays x 256 MiB: far beyond any cache
    // Peak = cores x GHz x FMA units x lanes x 2 flops. Fill in your machine.
    let (cores, ghz, fma_units, lanes) = (threads as f64, 3.5, 2.0, 8.0); // e.g. AVX2 f32
    let peak = cores * ghz * fma_units * lanes * 2.0;
    let ridge = peak / bw;
    println!("bandwidth {bw:.1} GB/s, peak {peak:.0} GFLOP/s, ridge {ridge:.1} flop/B");
    // (name, flops per element, DRAM bytes per element)
    let kernels = [
        ("saxpy y=a*x+y f32", 2.0, 12.0),
        ("dot f32", 2.0, 8.0),
        ("stencil 5-pt f32, perfect reuse", 5.0, 12.0), // read 4, write 4, write-allocate 4
        ("matmul f32, 128x128 tiles", 32.0, 1.0),       // AI = b/4 for tile size b
    ];
    for (name, flops, bytes) in kernels {
        let ai = flops / bytes;
        let p = attainable(peak, bw, ai);
        let bound = if ai < ridge { "memory" } else { "compute" };
        println!("{name:34} AI {ai:6.2} -> {p:7.1} GFLOP/s ({:4.1}% of peak, {bound}-bound)", 100.0 * p / peak);
    }
}

On an M4 Pro this measures about 144 GB/s with triad accounting. Apple's figure is 273 GB/s.

Complexity and constants

• Bytes per element for common kernels (f32): copy 8 (12 with write-allocate), axpy 12, dot 8, sum 4, elementwise c = a op b 12-16.
• Hierarchical rooflines: L1 and L2 bandwidth are order of 5-20× DRAM per core. A kernel that is DRAM-bound may become compute-bound once tiled for L2.
• GPU: H100 SXM HBM3 3.35 TB/s, 989 TFLOP/s dense BF16 tensor. Ridge ≈ 295 flop/B, so only large GEMMs and fused attention reach the compute roof.

Efficiency

• Cost: the analysis is O(1): count flops and DRAM bytes, take min(peak, AI × bandwidth). On the page's example machine (16 cores AVX2, dual-channel DDR5-4800: 1792 GFLOP/s, 76.8 GB/s, ridge 23 flop/B) the bounds are 12.8 GFLOP/s for axpy, 38 for untiled matmul, 1792 for 128×128 tiled matmul. H100 SXM BF16: 989 TFLOP/s over 3.35 TB/s, ridge ≈ 295.
• Levers:
  • tiling for cache — AI 0.5 → 32 (b/4, b = 128, example above) — kernels with data reuse such as matmul
  • kernel fusion — fewer DRAM round trips — chains of elementwise ops, attention
  • lower precision (BF16, int8) — fewer bytes per flop — memory-bound kernels
  • batching — AI ≈ 1 → ≈ B for LLM decode — matrix-vector inference
• Naive → standard → state of the art:
• naive: untiled i-k-j matmul, n ≫ cache · 38 GFLOP/s bound, 2 % of peak
• standard: tiled 128×128 matmul · 1792 GFLOP/s, at the compute roof
• state of the art: fused, tiled, low-precision GPU kernels (FlashAttention, CUTLASS) · reach the 989 TFLOP/s BF16 roof only for large GEMMs and fused attention

Implementation notes

• Raise AI by: tiling for cache reuse, fusing consecutive elementwise ops into one pass (kernel fusion, flash attention), smaller types (f32→bf16 halves bytes; int8 quarters them), compressing data, or computing values instead of loading them.
• Raise the memory roof by using all cores, NUMA-local memory, prefetch-friendly access, and non-temporal stores for write-only outputs.
• Raise the compute roof with SIMD and enough independent accumulators to cover FMA latency.
• Measure, don't guess: Intel Advisor and NVIDIA Nsight Compute draw rooflines from hardware counters. LIKWID (likwid-perfctr) reports flops and bytes on Linux.

Pitfalls

• Counting algorithmic bytes and forgetting write-allocate and cache-line granularity (touching 4 B of each 64 B line moves 64 B).
• Using the datasheet bandwidth. Achievable is lower, and one core alone gets only a fraction of it.
• Peak from the wrong instruction mix. Only FMA-saturated code reaches 2 flops per lane per cycle. Adds or multiplies alone give half.
• Applying the DRAM roof to data that fits in cache: use the cache-level roof instead.

Variants and alternatives

In production

Intel Advisor, NVIDIA Nsight Compute, LIKWID, and NERSC's roofline methodology for HPC codes. GPU kernel design (FlashAttention, CUTLASS tile sizes) is driven by this arithmetic.

Implementations

Sources

• Samuel Williams, Andrew Waterman, David Patterson, Roofline: An Insightful Visual Performance Model for Multicore Architectures, Communications of the ACM, 2009. paper
• Aleksandar Ilic, Frederico Pratas, Leonel Sousa, Cache-aware Roofline model: Upgrading the loft, IEEE Computer Architecture Letters, 2014.
• John D. McCalpin, Memory Bandwidth and Machine Balance in Current High Performance Computers, IEEE TCCA Newsletter, 1995. paper
• John D. McCalpin, STREAM: Sustained Memory Bandwidth in High Performance Computers. site
• NVIDIA, H100 Tensor Core GPU datasheet. page